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A NOTE ON UNIQUENESS PROOFS FOR BOUNDARY- 
VALUE PROBLEMS IN POTENTIAL THEORY AND 
STEADY HEAT CONDUCTION 
By M. E. RAYNER (Westfield College, London) 

Received 22 May 1952] 


SUMMARY 
For boundary-value problems in an infinite region in two-dimensional potential 
theory and steady heat conduction, the logarithmic singularities of solutions, in 
general, prevent our proving uniqueness. If, however, we allow the boundary to 
tend to infinity in two directions, and consider the problem in a semi-infinite region, 
for a very wide class of boundary conditions, we can find additional information 
which enables us to prove uniqueness (except possibly for an arbitrary constant) 


by classical methods. 


1. Introduction 
THE classical uniqueness proof for boundary-value problems in potential 
theory is based on the Green identity (1). The argument can be applied 
without difficulty to finite regions, and it can be extended to infinite regions 
in space provided that the potential is regular at infinity.t This condition 
is satisfied for fields due to point sources, but it fails to be satisfied for fields 
with plane symmetry due to line sources, which give rise to problems of 
‘two-dimensional’ potential theory. For these problems the elementary 
potential is of the form A logr, and thus tends to infinity as r tends to 
infinity, while in the three-dimensional problem the elementary potential 
is of the form A/r, and thus tends to zero as r tends to infinity. In either 
case, the constant A is determined by the source strength. 

Ifin the two-dimensional case, at large distances, we subject the potential 
v to be of the form Alogr+f(0)+9(r, 4) 
for some undetermined constant A and functions f(#) and g(r,@), where, 
however, g(r, @) is regular at infinity, and prescribe either 


(i) its value, 


3 4 C 
(ii) the value of | —», 
on 


(iii) the value of : vthey (h>0) 
on 


+ A function g(r, , where all the 6’s are angles, is said to be regular at infinity if 
g(r, 9,...) = O(1/r), and (é/ér)q(r, 8,...) = O(1/r?) as roo. 

{ We adopt the convention that @/én denctes differentiation along the outward normal 
to any region. 


[Quart. Journ Mech. and Applied Math., Vol. VI, Pt. 4 (1953)] 
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on a closed curve C, this being the interior boundary of an infinite region, 
then in cases (i) and (iii) v is not always uniquely determined. In fact. 
suppose that C is the circle r = 1, and that v*, supposed to be the difference 
between two solutions with prescribed values on C, is given by v* 





log r. 
Then v* = 0 on C, and satisfies the condition at infinity, but it does not 
vanish identically. Hence the first boundary-value problem is not unique, 
Again, suppose that C is the circle r = a, where a (> 1) is the unique 
solution of the equation 


aloga L/h, 
and let v* = logr. Then 
O os * 
—v*+ hv’ —l/r+hlogr 
on 


vanishes for r = a, but it does not vanish identically, and hence the third 
boundary-value problem is not unique. 

In case (ii) the proof of the uniqueness of v requires an additional argu- 
ment involving Gauss’s theorem. It should be noted that the appeal to 
this theorem is not necessary in the corresponding three-dimensional 
uniqueness proof. 

In a similar manner it may be shown that the boundary-value problem 
for composite regions separated by a closed curve, such as the case of two 
dielectrics, or that of steady heat flow in two media with or without contact 
resistance at the boundary, have unique solutions except for arbitrary 
constants. 

This note is primarily concerned with the uniqueness proof for semi- 
infinite regions. The question seems to deserve some attention for two 
reasons. First, despite what has been said above, uniqueness (except 
possibly for an arbitrary constant) can be proved for semi-infinite regions 
for all the usual boundary conditions. Secondly, in the literature on 
potential theory, a uniqueness proof for semi-infinite regions is usually 
omitted,} although one of the simplest examples for which the boundary- 
value problem can be solved analytically is the half-plane bounded by a 
straight line. I believe that the present note therefore fills a gap. 

In addition, the proof provides information about the total source 
strength of the so-called ‘image systems’ for general semi-infinite regions. 
This information is also of interest in the three-dimensional problem, 
although in this case it is not required for the uniqueness proof. 

We first treat the very simple cases of the first and second boundary-value 
problems, which are closely related to the image methods of electrical con- 
ductor problems and hydrodynamical problems of irrotational flow of an 


+t But see, however, H. and B. 8. Jeffreys, Methods of Mathematical Physics (Cambridge, 
1950), section 6.074. 
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inviscid fluid. Afterwards we deal with the rather general problem of steady 
heat flow in a composite region with contact resistance. The proof of the 
uniqueness of the solution of the third boundary-value problem or the 


dielectric case can be constructed on similar lines. 


2. Notation 


We shall suppose throughout that the whole plane is divided into two 
semi-infinite regions G, and G, by a regular curve I’ which goes to infinity 
in two directions.} If O is any point in G,, and P, and P, are points on IT, 
let w be the angle P, OP,. As P, and P, go to infinity in different directions, 


suppose the limit of the angle w is w,. Then w, = 27—wy. 


3. First boundary-value problem 

Suppose v’ and »” are two solutions of Laplace’s equation in G,, which 
assume a given value on’. We will suppose that, at large distances, v’ and 
r’ are of the form A’logr+f’(@)+g'(r,@) and A” logr+f"(@)+ 9’(r, 8). 
respectively, where g’(r,@) and g”(r,@) are regular at infinity. Then 


y'—y" 0 on [., so that as r >oo on T, 
(A’— A" log r-+-(f’(6)—f"(8)) = O(1/r). 


Therefore A’ = A”. Thus any two solutions can differ only by a term 
which is O(1) as r > 00, and by the classical uniqueness theorem it can be 
shown that these solutions can, in fact, differ only by a constant; but since 
the solutions take prescribed values on T’, the solutions are identical, and 
so the first boundary-value problem has a unique solution. 

We can look at this problem in a slightly different way. Suppose that 
there is an electric field with potential v, due to charges which are all in a 
finite sub-region of G,. If we now introduce a conducting, infinitely long 
sheet along I’, how will the field be altered? Ifthe perturbing potential is w, 
then the modified potential is v,+u and so u v, on I’, from which it 
follows that for large values of r, u —A,logr+O(1)ifv, = A,logr+O(1) 
asr—> oo. If we now interpret the perturbing potential wu as due to a system 
of image sources situated in G,, we have the result that the algebraic sum 
of these sources is equal in magnitude but opposite in sign to the sum of 
the sources of Uo- 

It is evident, by the same argument, that this result also holds in three 
dimensions, where all the sources of vy lie to one side of an infinite open 


regular surface. 


+ Thus exclud for instance, the case of the curve formed by the two spirals r 0 


and G x, Wwuere VU x ons 
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4. Second boundary-value problem 

Suppose that the normal derivative of the potential is prescribed on [. 
and that v’ and v” are two solutions of Laplace’s equation in G,, which at 
large distances behave like 

A'logr+f'(@)+q'(r,6) and A” logr+f"’(@)+ g’(r, @) 
respectively. It is assumed that in G,, g’(r,@) and g"(r,0) are regular at 
infinity. 

By Gauss’s theorem, 

” c , ” 

v") ds (v' —v ) ds 0. (1) 

on : on 

TR) Si(R) 
where S,(R) is the are in G, of a circle S(R), radius R, and T'(R) is the 
finite are cut by S(R). 
Then as R + o, since (é/én)(v’ —v") 0 on [’, from (1), 
w,(A’—A”) > 0, 

i.e. A w= A. 
With this information the solution can be shown to be unique except for 
a constant. 

Also, we may consider the hydrodynamical problem when 7v, is the 
velocity potential of a flow due to a given system of sources in G,. To find 
the effect of introducing a wall along I’, we define the perturbing potential 
u, so that 
(Up+U) 0 


on I’. Suppose that vy behaves like A,logr+f,(@)+g)(r, 0), and w like 
Blog r+-f,(9)+9,(7, @), at large distances, where g,(r,@) and g,(7,@) are 
regular at infinity in G,; u can be interpreted as the velocity potential due 
to a system of image sources in G5. 

If we assume that all the sources in G, lie in some finite region, it is 
possible to find an R such that, by Gauss’s theorem, 


fad C . oO 

— (vy +u) ds = — (U9 +u) ds = 27A,. (2 

J on J on 

T(R)+Si(R) SiC?) 
" : r we 
Since lim (vp tu) ds = w,(A,+ B) 
4 0 1 ol ’ 
R-o , on 
S,(R) 
») 
; (27—w Ws : 
from (2) B= — 1) 4, = —* A). (3) 
Wy W, 


Thus the ratio of the algebraic sum of the image sources to the sum of the 
sources in G, is wo/w,. 
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A completely analogous argument in three dimensions shows that the 


algebraic sum of the source systems are related by (3), but in this case w, 


ind w, must be interpreted in terms of solid angles. 


5. Steady heat conduction in a composite region with contact 


resistance 


Let v, and v, be the temperature distributions in G, and G, respectively, 
then v, and v, are solutions of Laplace’s equation except at those points 


it which heat sources are situated; on T° they satisfy the conditions: 


C\ 1 H(v, v1). (4) 
C,—? = H(v,—1,), (5) 


where C,, C,, and H are positive constants, and ¢/cv denotes differentiation 
ulong the normal to [' from G, to G5. 


is the temperature distribution due to a given distribution 


Suppose Uo 
of sources which lie within some finite region of G,. To find the effect of 
introducing the boundary I’, we define perturbing temperature distributions 


uv, and u. so that 
i = Vo Un Uo. 


u, and uw, will have no singularities in G, and G, respectively. 

We assume that, for large r, vy has the form A, logr+f,(@)+g,(r. 4), 
u, the form Blog r+ f,(@)+9,(7, @), and u, the form Blog r+ f,(@)+g.(r.9), 
where g,(7, @) and g,(r, @) are regular at infinity in G, and g,(r, @) and g,(r, 4) 
are regular at infinity in G,. 

[t is possible to find an # so that S(#) contains all the singularities of vp, 
whence from Gauss’s theorem, 


é . 
(v,+u,) ds | (vp+-u,) ds 277A. (6) 
on 


« . 


S,(R rh) 


In the same way, if S,(2) is the are of S(R) lying in G,, 


. . 


c - 
(v,+us) ds | (vyt+us) ds = 0. (7) 
J on ‘i 
So(R Ti?) 


Using (4) and (5). from (6) and (7) it follows that 
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When R — o, the left-hand side of (8) tends to 
C, (Ag+ By) +C, w(Ao+ By), 
while the right-hand side is constant. Therefore 
(C,—C,)w, Ap +C, w, Bi +C,w, B, = 0. (9) 
teverting to (7), and using (5), 


a t . 


é H 
| &m (vp tus) ds+ G. 


2¢ 
S2(R) r(R) 


(u,;—Uy,) ds = 0. (10) 


In the limit as R — 00, the first term on the left-hand side of (10) converges 
to (A,+ B,)ws, therefore the second term converges. This implies that as 
R + «, along I, U4 — Uy -> 0; 
thus B, = &,; (11) 
and from (9) and (11), we have 

w,(C,—C,) 


1 Oy + a2, 


B, = B, = Ay. (12) 
From (12) it follows that for the given conditions on an infinite boun- 

dary, and a given distribution of singularities, the perturbing temperature 

distributions can be proved to be unique except for an arbitrary constant. 
For a straight-line boundary, w, = w, = 7, and so 


_(G—G) 


B, = B,= C4 C, As. 


In three dimensions, when v, has the form A,/r+-g(r, 0,4) where 
g(r, 8, d) O(1/r?) and - g(r,9,¢) = O(1/r°) 
cr 


as roo, and with corresponding forms for u, and u,, it can be shown, 
by an analogous argument, that (12) remains true when w, and w, are inter- 
preted in terms of solid angles. 


REFERENCE 
1. O. D. KELLoGG, Foundations of Potential Theory (Berlin, 1929), 211—15. 
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NOTE ON LIN’S ITERATION PROCESS FOR THE 
EXTRACTION OF COMPLEX ROOTS OF ALGEBRAIC 
EQUATIONS 
By J. MORRISt and J. W. HEADY 
Received 8 May 1952] 


SUMMARY 
\ simple iterative procedure for the extraction of factors of algebraic polynomials 


has been given by Lin (1). The rational basis of this procedure and the conditions 


for its convergence are here discussed. 


Introduction 

AMONG many possible methods of determining quadratic factors of algebraic 
polynomials, that due to Lin (1) has been found particularly useful by one 
of the authors in the solution of practical problems connected with servo- 
mechanisms. The weakness of Lin’s method is the lack of specific criteria 
for its convergence. 

The rational basis of Lin’s process is here investigated, and clear-cut 
criteria for its convergence are given. These criteria are for the most part 
qualitative, because the appropriate analytical expressions contain the 
unknowns’ being sought and the criteria found are only truly rigorous in 
the immediate vicinity of the roots. 

The extraction of a quadratic factor from a quartic equation and the 
simplest case of extracting a linear factor are discussed in detail; the general 
case of extracting a factor of degree r from an equation of degree n is briefly 


considered. 


1. Description of the process 
As an illustration of the procedure consider the quartic 
F(x) xt+ 23+ 8x?+ 10a-+-7. 
In the absence of any indication of an appropriate approximation to a 
quadratic factor, we take the last three terms of the quartic, viz. 
8a?+10x+7, 
divide this coarsely by the coefficient of x? and thus make 2?+-1-25x-+-1, 
for example, a starting factor. Then we proceed as follows: 
xv? 1-25a-+-1)at+ 223+ 827?4+ 102+ 7(x?+0-75x 
4) 1-25a3 | a 
+ 7x24 10x 
0-75a3 + 0-94a?2+-0°75a 
6-06a2-+ 9-254a-+7 
Research Consultant. t B.B.C. Research Department. 


(Quart. Journ. Mech. and Applied Math., Vol. VI, Pt. 4 (1953)] 
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HEAD 


and thereby obtain the quadratic remainder 





6-06x7+-9-25a+7 = 6-06(a?+ 1-5382+- 1-155). 

































Then we take the factor x?+-1-538x-+-1-155 for the next approximation 
and repeat the process with the original quartic; and so on.+ Thus we 
obtain in succession 
F(a) = (a2-+1-538x-+ 1-155) (22+ 0-462) + 6-129(a2-+ 1-544a-+- 1-142), 
F(x) = (v?+-1-544a+ 1-142)(a?+ 0-4562) + 6-154(a?+ 1-5402-+- 1-137), 
F(x) (xv? +-1-540x-+- 1-137)(a?+-0-4602) + 6-156(a?+ 1-540x+- 1-137), 
Hence after four iterations we derive the quadratic factor 
xv? + 1-540a2-+- 1-137 
and, in the course of the procedure, the remaining factor 
x?-+0-460a2-+ 6-156, 
both computed to four figures. As a check we notice that 
(v?-+- 1-540a-+- 1-137)(2?+- 0-460a2-+ 6-156) 
x*+-223+-8-0012?-+- 9-998a-+ 6-999. 
It may be remarked that in the early stages, and until the factor we are 
trying to extract appears to have become stabilized, it is unnecessary to 
work to any but superficial ‘accuracy’ for the reason that in a fully iterative 
process errors are only of transitory effect until the final stages, when we 
must have regard for arithmetical accuracy to the required number of 
significant figures. 


2. Convergence of the process 

We notice in the numerical example given in the preceding section that, 
for the particular quartic considered, the successive approximation process 
for the extraction of a quadratic factor was operated rapidly, but such 
quick convergence will not necessarily occur in other cases. As in all 
iterative processes, criteria specifying the conditions for convergence 
become of paramount importance if the process is to be of practical value 
in actual arithmetical computation. 

Reverting to the particular numerical example worked out, we observe 
that for a root x of the equation F(x) 0 


(v?+ 1-25a-+-1)(a?+- 0-752) — 6:-06(a?-+ 1-538a-+ 1-155), 
(~?-- 1-538a-+ 1-155)(x?+ 0-4622) — 6§-129(27?+- 1-5447-++- 1-142), 
(v?+-1-544a-+ 1:142)(7?+ 0-45627) -6§-154(a?-+- 1-540a-+ 1-137). 


In Lin’s original paper the process by which the successive divisors are obtained is 
described somewhat differently, but the successive divisors themselves are identical with 
those given here. 








Thus 


ratio: 


and 
ratic 
bein 
It 
fact 
for 1 
thai 
abo 


T 


and 


Eq 


Le 








ition 


Ss we 














NOTE 





ON LIN’S ITERATION PROCESS 


Hence 


y2+ 1-540a-+- 1-137 
— 9 
W402) tz 


04622) (a?- O-606e) » 
4 


—_— 1-125¢7-++-1). 
6-129 6-154 


6-06 


Thus the convergence of the process in this case depends on the successive 


ratios 
| 
r(a 0-75) u(x- 0-462 ux 0-456) 
) ) a) 
R, oe, ———, base, 
6-06 6-129 6-154 


and so on, becoming numerically less than unity; and the smaller these 
ratios the more rapid the convergence towards the specific quadratic factor 


being sought. 





It seems therefore that starting with a more or less arbitrary quadratic 
factor there will eventually emerge from the process that factor, if any, 
for whose roots the foregoing ratios R; tend to become numerically smaller 
than unity. Otherwise the process will be divergent and in consequence 
ibortive so far as producing a quadratic factor is concerned. 


To rationalize the process for the quartic in general let 
Py j 3 2 . ” 
F(x) ? A, 2° +-A,2°+a4,%+Qy; (1) 


ind let this be written 





| a, 
F(x) (2*-+ p,a 7,.)(2x* P. x) >. (x? +-p, 12 +p) (2) 
Grit 
Equating corresponding coefficients we must have 
p F As (3) 
3 ] 
p, P+, = Ag—Ay (4) 
r+ 
) ~ 
Y i ay ie (9) 
Gr4t 


Let the quadratic factor to emerge from the application of the process 


assumed convergent for that factor) be a?+ pa-+-q, so that 


’ l > . 
F(a) (a pa q)(a* Px) +. (a2 pa q) (6) 
q 
Also let us write 
P P-+-Sp; Uy GT Dr od = P+ X,, 


where €.. X,, are assumed to be small quantities of the first order. 


(5), and neglecting quantities of higher order 


(7) 


Substituting these in (3), (4), 


than the first, we find ; 


ce XxX. 0, 
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Ao 


é.P+X,p- Ny = 5 Urey (8) 
q 
‘ a Ag Pp 
X,q+7,P = — Gert a panes (9) 


[t must again be pointed out that (7), (8), and (9) are only rigorously true 
in the immediate neighbourhood of the roots, and that this neighbourhood 
is unknown initially. 

In view of the relation P = a,—p we obtain from the relations (7), (8), 
(9) the two equations in €,, 7, as follows: 


a 
[plas 2p)-+ q\f,—(4gs—2p)n, = qe (10) 
(a,—2p)é, +», ? Nit (11) 
Now let &. = &N’, Nr = 7’," (12) 


where A is a parameter and &, 7 are appropriate constants. 

Substituting these expressions for €,, », in the equations (10) and (11), we 
derive two latent root equations in A as follows: 
E—(a,— 2p) = 0, (13) 


[Pla —2p)- q—“r 
q 





~ a 
¢ 0 4 
(ay 2p)8-+(1 -- 2A) = Q. (14) 
q* 
The corresponding characteristic equation in A is thus the quadratic 
2 9 
4 2__ %o P41 p(a3— 2p) 
g 4 q 
and for the convergence of the iteration process under consideration as 
applied to the quartic (1) both roots of (15) in A must have modulus less 
than unity.+ For the quartic of section 1, we have ag = 2, a, = 7, and for 
the factor x?+-1-540x7+ 1-137, we have also p = 1-540; q = 1-137. In this 
case (15) reduces to 


js [(as—p)(ag—2p)+q] = 9; (15) 


33°34A2— 3-307A+. 0-6402 = 0. (16) 
The roots of this quadratic are complex and have modulus 
(0-6402/33-34)! = 0-139. 

This gives the rational explanation for the rapidity of the convergence 
of the process in the particular case considered. 

For the other factor, viz. x?+-0-460x+- 6-156, we find that the modulus 
of the roots is 6-10 and thus this factor could not be found primarily by the 
process in the ordinary way. 

+ If we eliminate x between the equations x?+ px+q = 0 and x?+ Px-+ (Aay/q) = 0, we 
obtain (15). It follows that if 2, and 2, are the roots of z?+-pa+q = 0, the roots of (15) are 


-q(xj+-Px,)/ay and —q(x3+ Px,)/ay. The significance of the ratios called R; in section | is 
thus explained. 
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We notice that the solution of the iteration equations (13), (14) above is 


E, = EN+E243, (17) 


(9 Vr AL+ AeA, (18) 


— where A,, A, are the twe roots in A of the quadratic equation (15) and 
‘ 2 


ail Ei: €, 2 are appropriate constants. 

Lin (1) suggested that when his process was not convergent the given 
| equation should be transformed, either by squaring the roots or by increas- 
| ing or decreasing them all equally. He gave examples for which such 
transformations did produce convergence. In other cases, however, the 
(10 convergence may be retarded; for example if the roots of the equation 
discussed in section 1 are squared, the equation becomes 

y*+ 12y'+ 38y?+ 12y+49 = 0 (19) 


12 ° 
\ and the convergence is retarded. 





If as a result of any transformation (1) becomes 
), We , , 

G(y) y* +b, y+, y? +5, y+-5,, (20) 
‘n approximate equation corresponding to (15) can be formed by replacing 
in (15) 

ay by by; as by 63; p by (b,/b,);  ¢ by (b,/b,). 

(14) Unless the roots of this approximate equation have moduli appreciably less 
than unity, the right transformation has not been found. If better approxi- 
mations to p, g than (6,/b,) and (b,/b,) are known, they should clearly be 


(15 used. 





n as 3. The general case 
less 
d for 
this 


Although our discussion here is limited to the extraction of a quadratic 
factor from a quartic, Lin’s method can be used (when convergent) to obtain 
a factor of degree r from an equation of degree n (r < n). If the n-ic is 


(16 atta, _,x"1+...+a,22+a,r+a, = 0 (21) 
our first divisor is, in effect, 


‘ r , 99 
(z/a,)(a,a"+a,_, 2"-!+-...+a,2%+4,); (22) 
ence ‘o . 7 . . . . . 
and if the remainder after this division is 


[ulus B(x) = b,a7+6,_,277-1+...+6,24+ bo, (23) 


y the the second divisor is x B(x)/b., and so on. 


Criteria for convergence are obtained as for the quartic, but the equation 
0, we corresponding to (15) becomes increasingly complicated as r and n increase, 


ner and for convergence all its roots, whether real or complex, must have 
n 18 
modulus less than 1. 
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4. The special case of a linear factor 
The extraction of a linear factor has interesting special features. (on- 
sider, for example, the simple cubic 
F(x) = 23+ 622+ llz--6. 24) 


We obtain in succession 


F(x) = (w+ 2-67)(a?+-3-332)+-2-11(a-+ 2-84), 

F(a) (~-++ 2-84)(x?+- 3-162) + 2-03(a-+ 2-96), 

F(x) (x+ 2-96)(a?+-3-042)+ 2-00(a-+-3). (25) 
Now F(x) = (w#+1)(a+2)(a+3), (26) 


and it will be found that the process is convergent for the factors (a+1) and 
(x+-3) but not for the factor (a-+-2). 


In general let F(x) represent a polynomial of any degree and let 


F(x) = (x—p,)xP,(x) “0 (x—p,44)- (27) 
Pr+i 
Then F(p,) = - “0 (P-—Pr+1)s (28) 
Pr+i 
or 4 _ WF (P,) (29) 
Prot Pr 


which is the symbolical representation of Lin’s iterative process applied 
to the case of a linear factor. 
For the derivation of a convergence criterion let p, = p+-€,, where p is 


a root of F(x) = 0, and € is assumed to be a small quantity of the first order. 





' 5 ‘ »F’(p)]. : 
We thus find that Ci | es ce (30) 
A 
the corresponding A equation being 
PF’ 
ee Pe dl! od (31) 
Xo 
where &, = &r (32) 


is the solution of the ‘difference’ equation (7). 

Clearly the convergence will depend on the numerical value of A being 
less than unity and the €, will comprise a pure geometrical progression. 
In particular, if "a , = 33 

I F(x) (x-+ p,)(@+ pg)...(w+p,,) (do) 


for the factor («+-p,), the appropriate value of A obtained from (31) becomes 


Np,) = 1—\P2— PMPs Pr) Pn—Pr) | (34) 
P2P3---Pn 
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; If therefore F(x) is transformed by increasing each root by a, A(p,) will 
Con have to be replaced by 
(Po—P)( p: ),)...(9,— PD) _ 
: u(p) P2—-PiK\P3s—P\ } P} (35) 
(24 (Pot a)(P3+a)...(p, +a) 
and as « varies, w can vary over a wide range of values. 
5. Conclusions 


We have obtained definite convergence criteria for Lin’s process when 


a quadratic factor is extracted from a quartic, and the way in which a root 





(25 isreached when a linear factor is extracted by the process has been examined 
sa for an equation of general degree whose factors are known. If the equation 
- is transformed, the criteria can be applied to the transformed equation, 
pe but any particular transformation, such as squaring or increasing the roots, 
, can either improve or retard the convergence. The effect of any trans- 
formation on convergence in a particular case can be roughly gauged by 
= applying the convergence criteria to the transformed equation, with the 
. unknown roots involved in those convergence criteria replaced by the best 
available approximations to them. 
s 
, REFERENCE 
29 1, SutH-NGE Lin, ‘A method of successive approximations of evaluating the real 
and complex roots of cubic and higher order equations’, J. Math. and Phys. 20 
— 1941), 231-42, 
p is 
det 
i 
(30 | 
| 
| 
} 
3] 
(32 
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THEORETICAL CONSIDERATIONS ON FREE 
CONVECTION IN TUBES 
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SUMMARY 


Methods are developed for predicting the flow and heat transfer due to free con- 
vection in heated vertical tubes, closed at the bottom and opening into a reservoir 
of cool fluid at the top. These include methods of predicting whether the flow is 
laminar or turbulent, and whether the boundary layer of heated fluid fills or does not 
fill the tube, or fills it with a stagnant region near the closed end; and six methods 
of analysis are needed for the various combinations of these cases. The flow found 
depends on a certain modified Grashof number G, (see List of symbols) and on the 
length—radius ratio l/a. For laminar flow, indeed, it depends only on the quotient of 
these quantities (Fig. 7), but laminar flow becomes hard to achieve for large G,. For 
large entry disturbances turbulence is expected for G, greater than about 10‘, and 
then the shape of the heat-transfer curve changes with1/a. The theory of the turbulent 
flow, which is based on an assumption relating the exchange coefficient to its value 
for ordinary pipe flow, indicates that turbulence decreases heat transfer when the 
boundary layer fills the tube and increases it otherwise. This leads to a peculiar 
‘bucket’ in the heat-transfer curves (Fig. 9). 

The Grashof number at which the heat transfer becomes most effective rises 
steeply with l/a (Fig. 10). 

The theory is extended to apply to a tube closed at both ends, part heated and 
part cooled. There is some discussion of a mechanical engineering application in 
which the external force is not gravity but centrifugal. The whole theory is tentative 
only, and a critical evaluation is postponed until a larger body of experimenta! 
information is available, but some reference to experiment is included in a ‘Note 
added in proof’ at the end of the paper. 
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FREE CONVECTION IN TUBES 


List of symbols 


acceleration due to gravity 

externally produced acceleration (may or may not be g) 
coefficient of cubic expansion 

temperature 

temperature at walls 

temperature on the axis at the orifice 

length of tube 

radius of tube 

kinematic viscosity = p/p 

thermal diffusivity = k/pc,, 

density 

specific heat at constant pressure 

viscosity 

thermal conductivity 

modified Grashof number (or ‘Rayleigh number’) based on 
x(7,—7,) fF 


VK 


length 


«(Ty —T;) fa* 


modified Grashof number based on radius 


VK 
Prandtl number = v/« = pc,/k 
distance from closed end 
distance from axis of tube 
axial velocity (i.e. in X-direction) 
radial velocity (i.e. in R-direction) 
pressure 
non-dimensional axial velocity = a?U’/«l in laminar flow and 
KaU /v in turbulent flow 
non-dimensional radial velocity = aV /« (in laminar flow only) 
non-dimensional temperature = afa*(7,—T7)/v«kl in laminar 


flow and vfa?(T, — T)/v? in turbulent flow 
non-dimensional distance from axis Ria 
non-dimensional distance from closed end X/l for laminar 
flow and K?X /a for turbulent flow 
value of ¢ (laminar flow) on axis at orifice (a/l)G, = (a/l)'G, 
total rate of heat transfer from tube 
Nusselt number based on radius = Q/2zlk(T,—T,) 
Nusselt number based on length = Q/2zak(T—T’) 
a function of x (with different significances in sections 4, 5, and 6) 
functions of x (section 4) 
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q mean rate of inward radial heat flow per unit area per unit time, 
divided by pc,, 

Jo value of q at the wall 

T mean shearing stress retarding upward flow from inside, divided 
by p 

T% skin friction divided by p, i.e. minus the value of 7 at the wall 

Tos average temperature over a cross-section 

A logarithmic decrement of turbulent similarity solution 

€ exchange coefficient (eddy kinematic viscosity) 

wv’ turbulent shearing stress (divided by p) = «@U/eR 

K,A constants in empirical equation for mean velocity in turbulent 
pipe flow 

R Reynolds number based on skin friction = av7,/v 

Q non-dimensional form of local heat transfer rate fa"qo 

Kv? v7, 

B,C, D functions of x (sections 10, 11) 

Z., t, temperatures (actual and non-dimensional) on the axis 

Nx Nusselt number based on heat transfer from a length X of tube 


measured from the closed end, and on the length X, 


x 
{| a» ax] /i(2y— TL} 


B, value of B where turbulence disappears 
H suffix denoting heated portion of tube 
C suffix denoting cooled portion of tube 


1. Introduction 
THE subject of this paper is the effect of confining walls (and floors) on 
laminar or turbulent free convection flow. But the author’s interest in this 
general subject arose through consideration of a specific practical problem. 

One would hardly have expected heat transfer by free convection currents 
to be useful in mechanical engineering, largely because the appropriate 
scales of length and time are too small for the acceleration due to gravity 
to have much effect. But free convection can take place in any external 
force field, and in particular the centrifugal accelerations (which in rotary 
machinery may be of the order of 104g) can produce very rapid free con- 
vection currents, causing extremely effective heat transfer. 

E. Schmidt first proposed to use these for cooling turbine blades, and a 
turbine incorporating his ideas was constructed in Germany during the war. 


Each blade contained a number of thin cylindrical cavities pointing radi- 
ally outwards from a reservoir of cool fluid in the hub. The effect of the 
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centrifugal force field was to cause the fluid heated at the circumference of 


' each cavity to move towards the hub and be replaced by cool fluid moving 


outwards down the centre. Schmidt has recently written an account (1) 
of this and later work. 

The parameter of greatest significance in heat transfer by free convection 
is the modified Grashof numbery (sometimes called Rayleigh number) 
(T,—1)) 0 


VK 


(1) 


Here « is the coefficient of cubic expansion, 7,— 7, the temperature differ- 
ence, f the external acceleration (usually g, but in the present application 
centrifugal and many times gq), / a typical length, v the kinematic viscosity, 
ind « the thermal diffusivity. (We shall sometimes base the modified 
Grashof number on the length / of the tube and call it G,, and sometimes 
on the radius a and call it G,.) The Nusselt number 
Q 
iS(h—T) id’ 


where Q is the rate of heat transfer across area S, k is the thermal conduc- 


(2) 


tivity of the fluid, and dis a typical length, is found to increase as G increases. 
For example, for heat transfer from a vertical flat plate with laminar flow 
NV o G (Pohlhausen quoted in (4)) and the increase appears to be rather 
more rapid when turbulent flow occurs, namely for G greater than about 10° 
Saunders (2)). 

\n important feature of Schmidt’s proposals was to maximize G not only 

by using a large value of the acceleration f but also by adjusting the pressure 
nd temperature of the cooling fluid as closely as possible to the critical 
condition, in which the viscosity is a minimum and the coefficient of expan- 
sion is very large.{ By these means it is expected to obtain modified Grashof 
numbers G, (based on the length of the cavities) of the order of 10%. 

A serious difficulty in the original system proposed by Schmidt is that 
sediment may accumulate at the ends of the cavities as a result of centrifugal 
separation. This is hard to avoid because the fluid in the cavities forms part 
of a large volume of fluid circulating between the turbine and a heat 
exchanger outside it. Therefore an alternative arrangement has been pro- 
posed (Cohen (1)) in which the cavities are completely enclosed; they 
stretch into, and are immersed in, a secondary circulating fluid in the hub 


The ordinary Grashof number is the same but with v? instead of vx in the denominator. 
Saunders (2) and McAdams (3, p. 248) showed that the influence of v and «x on free con- 

m flow is about equally important and that dynamic similarity between flows is much 

accurately achieved if the modified Grashof numbers are identical than if the Grashof 

bers proper are identical. 

The specific heat is also great in this condition, causing k to be large although the 
thermometric conductivity « is small like v. 


5092.21 D d 
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(Fig. 1). These concealed cavities then act like streaks of very highly 
conducting material within the blades, transferring heat from the blades 
to the cooling fluid. 

As a result of the conclusion of hostilities Schmidt’s turbine was aban- 
doned before it had been successfully run, although preliminary experiments 
had indicated that it was a workable proposition. However, the problem 









































of constructing a practical gas turbine incor- 
Y porating the idea is still being actively tackled. 
(on —_— Of the many difficulties, one is simply the 

a one Cavity - . 
4 only shown) very limited extent of knowledge concerning 
the actual flow and heat transfer in the cavities, 
There are so many variables involved, and each 
may take such a wide range of values, that a 
VW few experimental results are of dubious worth 
SY UY in the absence of theoretical considerations by 
Ht 7 means of which they may be extrapolated 
4 through the range of possible conditions under 
HY ; which the system might be used. It is the object 
d | Hub of this paper to provide such theoretical con- 
H 4 siderations; unfortunately very few experi- 
; i} } mental checks are possible at this stage, but it 
; ——— is hoped that a critical study may later be 
4 possible, since the volume of experimental 
a information is rapidly growing. If borne out 


Fic. 1. Arrangement pro- in broad detail, the considerations may be of 

posed by Cohen. help towards guiding the development of the 
technique, by reducing one at least of the many uncertainties in the way 
of progress. 

[t is also hoped that they may be of interest to scientists as illuminating 
a rather unusual kind of fluid motion, whose essential difference from 
problems of free convection which have been studied theoretically is in 
the effect of the constraining walls. 

The case of laminar flow is treated below at considerable length, before 
that of turbulent flow. From a scientific point of view this is natural since 
both subjects are equally interesting. From an engineering point of view 
it is equally desirable, even though in the practical application the flow 
will always be turbulent. For study of the corresponding laminar flow, 
being more reliable, is essential to bring out certain broad qualitative 
details of what happens,+ and also to give lower limiting conditions into 

+ Thus it will be found that there are three geometrically different régimes of flow, 
depending on the range of values in which the ratio of tube length to tube radius lies. 
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which the turbulent flow must subside at the Grashof number of transition. 
Further, it is quite possible that the laminar flow may have practical 
importance in problems of ordinary free convection under gravity. 

As the basic problem for throwing light on the subject we consider fluid 
ina cylindrical cavity, closed at one end and opening into a reservoir of cool 
fluid at the other end, as in Schmidt’s original proposal. The external 
acceleration f acts axially towards the closed end, and the solid boundaries 
are maintained at a given temperature 7). The fluid entering the tube along 
the axis at the open end is supposed to have the temperature 7; of the 
reservoir. The laminar flow is studied by means of the equations of motion 
in boundary layers. Approximate methods of solution are used in which 
the equations for momentum and energy transfer across a cross-section play 
a part. Such a procedure, based on the Karman—Pohlhausen method 
(5, ch. 4) in ordinary boundary layer theory, was first applied to free 
convection problems by Squire (5, p. 641). The calculation of the turbulent 
case is made possible by assumptions concerning the eddy viscosity, based 
on recent work by Reichardt (6) and by Townsend (7, 8). 

The basic problem just described will take up almost the whole paper, 
but it may be remarked at once that knowledge concerning this flow can 
be applied directly to predict conditions in the more complicated flow of 
Fig. 1. Here the hot fluid at the circumference of the cavities in the blade 
moves towards the hub and comes into head-on collision with the cold fluid, 
also at the circumference (where it is in contact with the coolant), which is 
moving outwards. We may suppose that the two mix perfectly and are 
forced in towards the centre. The temperature of the mixture, say 7,, will 
depend on the temperatures 7),,, and 7), of the blades and the coolant and 
on the ratio of the quantities of fluid colliding. Also this ratio will be equal 
to the ratio of the quantity of mixed fluid which moves out down the blade 
(along the axis of the cavity) to that which moves towards the hub. The 
flow in the heated part of the tube, on these assumptions, is identical with 
that in the ‘basic problem’ described above, if / signifies now the length of 
the heated part. That in the cooled part is the same with 7%, replacing 
7) and a new value for/. Finally the temperature 7, must be determined 
from the condition that it could be produced by mixing in the ratio in which 
the hot and cold fluid do mix at the join (which will itself depend on 7;).+ 
For example, if the two lengths were equal and the constants of the 
flow (including the external acceleration f) were uniform, 7, would, by 


+ We shall see that to join the solutions up it is sufficient to apply the conditions that 
the rate of heat transfer across a cross-section, and also the central temperature, are con- 
tinuous at the join; all other discontinuities being taken care of by the mixing process 
described above. 
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symmetry, be equal to 4(7%,+7% ). A more realistic case is considered at 
the end of the paper, and the above considerations applied in detail to it. 

The author expresses his gratitude to Mr. H. Cohen for introducing him 
to the problem and for assistance with it at all stages of the work. 


2. Laminar flow: the three régimes 

The characteristic non-dimensional parameters of the problem are the 
length-radius ratio l/a of the cylindrical cavity, the modified Grashof 
number G, or G, (it will appear that G, plays a slightly more important 
role; in any case G, = (//a)?G,), and the Prandtl number o = v/« (the ratio 
of diffusivities of shear and heat). It is not contemplated that fluids of 
Prandtl number differing markedly from 1 be used in the practical problem; 
this restriction will be used to simplify the calculations in the turbulent case, 
but the laminar flow calculations will be valid for all values of o except the 
very small ones typical only of liquid metals. 

When //a is sufficiently small (for a given G,, say) the flow up the sides 
of the cavity may approximate to the free convection flow up a vertical 
plate. This is a flow of boundary layer type, which was calculated in the 
laminar case by Pohlhausen. The experiments of Schmidt and Beckmann 
(4), and also of Saunders (2), gave excellent agreement with the calculations. 
For flow inside a circular cylinder it is known that the curvature of the wall 
as such does not affect the form or validity of the boundary layer approxi- 
mation. On the other hand, the solution would take no account of the 
confined space in which the motion occurs, which makes it necessary that 
the flux upwards at the walls be balanced by an equal flux downwards at 
the centre. In fact the Pohlhausen solution supposes no vertical velocity 
outside the boundary layer, the fluid entrained into the boundary layer 
being assumed to approach it horizontally. Hence the solution will be valid 
only when the cross-sectional area of the boundary layer is small compared 
with that of the whole cavity. This condition restricts greatly the length- 
radius ratio //a for a given G,. 

However, the solution in question is extended in section 4 by permitting 
a vertical flow outside the boundary layer in the centre of the cavity, as 
illustrated in Fig. 2. This yields the first of the three flow régimes. One 
might at first expect it to be valid up to a value of l/a for which the boundary 
layer first fills the whole tube. Actually, however, it should break down 
for a rather smaller value of //a, roughly that at which there is no longer 
a maximum in the volume flow of unheated fluid (i.e. of fluid outside the 
boundary layer) at the orifice cross-section. For it is only possible that 
unheated fluid be drawn into the boundary layer to become partially heated 
fluid, not vice versa. 
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When //a exceeds this critical value the boundary layer mixes with the 


central flow and, when a steady state is reached, fills the whole tube. The 
nature of the resulting flow is a little difficult to perceive intuitively, but 
an extreme case, with much larger l/a, can be visualized. This is one in 
which all tendency for the boundary layer to thicken with distance from 
the closed end has disappeared. Then the distributions of velocity and 











J : 
Fic. 2. Flow wit} Fic. 3. Similarity ré- 
boundary layer not gime with stagnant 


filling tube portion. 


temperature are similar at each section of the tube, only their scale increas- 
ing as the orifice is approached. It will be found that under these circum- 
stances the scales of both the temperature and (axial) velocity profiles 
increase linearly along the tube. The condition that the central temperature 
should rise from its value 7, at the orifice to the value 7, at the bottom 
then determines the length—radius ratio //a for this particular flow as a 
certain multiple of G, (section 5). 

This is the largest value of //a for which there is motion all along the tube. 
For greater values still (which we may regard as achieved by extending the 
tube beyond the closed end) the additional length is filled with fluid at rest 
at the temperature of the walls, the ‘similarity flow’ in the remainder being 
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unchanged. This third régime of flow, which occurs when //a exceeds the 
second critical value described above, is illustrated in Fig. 3. 

But a factor of order 10 separates the first critical value of l/a (corre- 
sponding to the breakdown of the boundary layer solution) from the second 
(corresponding to the onset of the similarity solution). The intermediate 
range must be filled by a solution which permits variation down the tube 
of velocity and temperature profiles which nevertheless fill the whole tube. 
Such a solution is obtained in section 6. The similarity profiles are attained 
near the orifice and more crimped profiles near the closed end, where the 
fluid rises in the middle as well as at the sides, falling only in an annulus 
between these two regions. 

The calculations for laminar flow are simplified by letting the Prandtl 
number o tend to infinity. Physically this corresponds to neglecting the 
inertia of the fluid compared with the forces due to its viscosity and the 
external acceleration f. It is known that this is a good approximation, even 
if c is really about 1, in the problem of free convection from a vertical plate, 
which is one extreme case of the present problem. The other extreme case 
(the similarity solution) is also worked out below for general o, and when 
o —> 0, and the difference found to be small when cis near 1. It is reasonable 
to conclude that the simplified calculations obtained in the limit as o > « 
are also a good approximation in all the intermediate cases, and, in these, 
only the simplified calculations are actually carried out below. 


3. The equations of the laminar flow 

The equations of motion of the problem are similar to the ordinary 
equations of free convection flow except that the pressure no longer takes 
its hydrostatic value. 

In all cases one is justified in simplifying the equations by the boundary 
layer approximation, in other words, neglecting the gradient of a quantity 
along the tube compared with its gradient along a radius. For in all cases 
either the flow is of boundary layer type (as in the extreme case of free 
convection from a vertical plate, for which the solutions of the boundary 
layer equations are known to be in excellent agreement with experiment) 
or else (when the boundary layer fills the whole tube) the length—radius 
ratio is large. 

With this approximation the equations of conservation of mass, heat 
and momentum respectively, in steady axisymmetrical flow, are 

aU ae V rT . wie «(ae aR) 


4 Vig pat, ae ee 3 
éX'eR'R ax? Vor = “lame Ror ?) 


OR? ' RéR eR 


ue = po aR ep a (4) 
4 p v 
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FREE CONVECTION IN TUBES 407 
Here X is measured along the tube from an origin at the closed end, and 
R radially outwards from the axis; the corresponding velocity components 
are UU and V. 7’ is the temperature, « the thermal diffusivity (thermometric 
conductivity), p the density, and p the pressure. 

In the momentum equation it is necessary to take into account the fact 
that p-! varies with temperature according to the equation 


P , Po | 1+ x( T’— T)|. (5) 


where suffix 0 signifies conditions at the wall R = a and «is the coefficient 


of cubic expansion. By this means the thermal buoyancy force is made 


apparent. For by (4), with R a (where U V 0), 
ap 227 au 
0 J ” +t val se : ’ (6) 
py dX OR? ROR] pa 


and hence, substituting for dp/dX from (6) in (4) and using (5), 
eu) al 
a 


eU aU 
eR?" ReR 


U (7) 


Equations (3) and (7) are three equations for U, V, and T' to be solved 
under the boundary conditions at the walls 
U V = Oand 7 = 7, for R 
and subject also to the condition that 7' 
X= 1, R = 0). 
The substitutions necessary to reduce simultaneously the equations of 


a and for X = 0, (8) 


T, at the orifice on the axis 


motion and the boundary conditions to non-dimensional forms, in which 
the minimum of parameters is retained, are 





i : ll 
. Kl , K os ae VK . 
l - 26, J v, 7 o——!> R=ar, X=Iz. (9) 
a? a xvfat 
Then equations (3) and (7) become 
cu ov v ct ot ort l ot 
4 { 0, u—-v- —+ ——, (10) 
On Or r Ox or or ror 
] ou ou o-u léu]" 
(uo) = t+ 4S, (11) 
Co Cx cl | or® ror}, 
while conditions (8) become 
uU=v=t 0 for r= 1 and forz= 0, (12) 


and the condition at the orifice becomes 


r=1,r=0 l 


4( 7 — 4 
it) vfa -T)) _ ; iL « (7) Gi. (13) 
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Thus (for a given co) the one parameter affecting the laminar flow régime 
is the parameter ¢, defined in (13) as the quotient 
modified Grashof number based on radius xvfa?(T,— T,) ve (14) 
length—radius ratio Lia ; 
This’ parameter determines which flow among the sequence described in 
section 2 actually occurs. 

Note also from equation (11) the simplification (section 2) achieved by 
allowing o to tend to infinity, and its physical significance as the neglect 
of inertia forces compared with buoyancy and viscous forces. 

Now the most important result of any calculation on the above equations 
will be the value of the total heat transfer rate, suitably put into non- 
dimensional form. For this purpose we use normally the Nusselt number 
based on wall area and tube radius, 


a 


l 
Q a for ‘ 
N, ———— 1X 15 
*  kKN—T7,)2xla/a (T%,—T,)l | ( 7). ”) 
0 
1 


1 / ’ 
| Z a dx, (16) 
ti. Cr}r=1 
0 


and occasionally that based on the tube length, namely N, = (//a)N,. The 
former, however, is a function of t, alone, and so the main results of the 
calculations will be graphs of N, against ¢,. But the distribution of heat 
transfer rate along the tube needs to be considered also. 

No attempt will be made below to solve equations (10) and (11) through- 
out the fluid, but only at the walls and on the axis and in the forms obtained 
by integrating them over each cross-sectional area, namely 


1 1 
; ar ot : 
| ru dr 0. | rut dr - . : (17) 
R Cx. cr r=1 
0 0 
1 1 
‘ee ae r 1 /ou cru 
| ru2 dr — | rt dr — — —— : (18) 
oOx | kj 2 \or 2, iat 


0 0 


These integrated forms express conservation of mass, heat and momentum 


equations (10) and (11) at the walls and on the axis are 


ot i] = | OF & re 
=+-—- 0, U et 
a Feri, OL] »=9 oT FOr}, 


] ou C7u sa). (19) 
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respectively for the fluid filling any cross-section of the tube. The forms of 
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Here the equations which would involve v are omitted, since by (17) and 
18) a solution for ¢ and wu alone is sufficient. The omitted equations result 
from the equation of continuity, which is fully taken account of by the 
first of equations (17). 

In sections 5 and 6 solutions of (17), (18) and (19) under the boundary 
conditions (12) and (13) will be obtained in the form of polynomials in r. 
Indeed, since by symmetry ¢ and w are evidently even functions of r, they 
will be taken as polynomials in r? (including terms up to 7*). This probably 
increases accuracy by assigning weight to different annuli according to 
their cross-sectional area. 

In section 4, where the ‘first régime’ of flow is studied, a solution with 
. variable boundary layer thickness is needed and (because of this added 
complication) a simpler form of the profiles is used, no attempt being made 
to apply the more delicate of the above conditions, namely those involving 
the values of second derivatives at the boundaries. In this we follow Squire’s 
work (5, p. 641) on free convection. For this purpose, an amended form of 
the momentum equation (18) is needed, in which (é?w/ér?),_, has been 
eliminated between it and (19). Since in this first régime of flow wu is inde- 
pendent of r outside the boundary layer there is no contribution to the 
square bracket in (19) from the upper limit, and the amended form of (18) is 


in ( onsequence 
1 
l ¢ 


1 
i 1/l Gu fou 
ruz di rt dr u - i : | ‘ (20) 
fom an ? 2Z\0 ox — OP J p= 


r 
0 0 


Except in calculating the similarity solution the terms in o~! in equations 


17) to (20) will be neglected, as explained in section 2. 


4, Flow with laminar boundary layer not filling the tube 

In this section the first régime of flow (section 2) is studied. In the extreme 
case of it when the boundary layer fills a negligible area of the tube, the 
flow is identical with the flow up a vertical plate. The full equations of this 
flow were solved exactly by Pohlhausen (4) with o = 0-733. He obtained 
N, = 0-517G} for the Nusselt number based on length / of plate, which gives 
the relation N, = 05174 (21) 
for the relation between Nusselt number based on tube radius and the 
parameter ¢, defined in (13). 

Squire (5, p. 641), on the other hand, gave an approximate procedure in 
which the velocity profile was taken as a cubic polynomial and the tempera- 
ture profile as a quadratic. His result, in our notation, is 


V 0-677(1 | iy (22) 
4 4d = 1° -- 


2lo 


\ j 
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For c = 0-733 this formula gives NV, = 0-550t}, which is only 6 per cent. in 
excess of the accurate value. 

We conclude that the error introduced by Squire’s approximate procedure 
is of tolerable order of magnitude. The additional! error due to replacing o-! 
by 0 depends, of course, on the value of c. In the above case it would be 
23 per cent., but this was for a low value of o typical of gases. For o = 2 
it would be only 10 per cent. These errors will have to be borne in mind 
in finally evaluating the theory (section 7). 

In this section the Squire procedure will be carried through with o-! = 0 
and with a vertical velocity permitted outside the boundary layer. Thus 
we assume that approximately 


f ty (0<r< B) 


= —_ 2 ? (24) 
r—B 5 Kee 
»|) = {1+6(r- " << 5) 
where f, y, 5 are functions of x. Equations (23) and (24) already satisfy the 
boundary conditions on r = 1. The equations (17) and (20) are now used 
to determine f, y, and 5. Thus the equation of conservation of mass in 
17) gives eye 2 
iaiaaias 5 — —>(3+2B +8) (25) 
3+ 2B8)(1—B)? . 
Equation (20) with o-! = 0 then gives 
ae (3+2 \(1 __ R\3 
oe so lo (26) 
36(3-+ 48+ 3B?) 
Finally the equation of conservation of heat in (17) gives 
d [(1—B)3(3- (45+ 1328+ 181B?+ 6283) Gi cn ] 
ae —ByPE tB) Bx pr. = t, —|[F(B)]| 
dx -30240(3-+.48-+ 3p?) 1—B 


da 
The complicated fraction F(f) is tabulated in Table 1, together with 


~~ 


— 
bo 
~~! 


B 
. Ff om te = 
4 | (1—f)dF(p)| . (28) 
The value of expression (28) at the orifice 2 = 1 is the parameter ¢,, so the 


table gives a relation between t, and the value of B at the orifice. The table 
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also gives values of the Nusselt number for given values of 8 at the orifice. 





This, by (15), is 
l 


1 
° 2t . 7 . 
| (Gag) @ = 26 FB) (29) 


l« 
0 


The cross-plot of NV, against t, obtained from Table 1 is shown in Fig. 7, 
where its asymptote NV, = 0-677¢ (the Squire solution for o-! = 0) is also 
shown. 

TABLE 1 











B B) t,/x (= t, at orifice) N, 
re) O01 455 1000 z715 
7 OO! 309 1290 3°37 
1094 16050 3°09 
3 000555 2400 4°12 
625 3770 471 
5 c 60740 5°50 
0 0002 35 14450 6°59 
7 ooII2 40500 g°O7 
5 0000 37 152000 13°5 
‘9 20000 5 2400000 27 
Io x x 





Although the solution exists formally for values of t, down to 1,060, the 
actual critical value of ¢, at which such a solution ceases to apply will be 
somewhat greater (see section 2), being determined by the condition that 
the flux of fluid outside the boundary layer ceases to be a maximum at the 
orifice. This flux by (24) and (26) is proportional to 





| 3+48+3p2 © 


6?(3+-B)(3+28)(1 —bP (30) 
which is a maximum for 8 = 0-38, corresponding to about t, = 3,400 and 
V = 4-6. This point is marked in Fig. 7 with a cross. If t, is significantly 
less than this value, then the solution requires that cold fluid should enter 
the boundary layer from outside and then leave it again without having 
been in the least affected. This is impossible and the solution must become 





increasingly inaccurate for t, < 3,400. 


5. The similarity régime (laminar flow) 

We pass now to the other extreme case and consider a solution of equa- 
tions (17), (18) and (19), under the boundary conditions (12), such that u 
and ¢ are proportional always to the same two functions of r, the factors of 
proportionality varying with x. A brief examination shows that the only 
possible form of variation with x is that uw and ¢ are both proportional to 2, 
measured from some origin. (For that they would have the same variation 
with 2 follows from (18), whence one concludes from either (17) or (19) that 
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it is linear. The conclusion is equally easily drawn from equations (10) 
and (11).) 

At the second critical value of t, = G,(a/l) described in section 2 this 
similarity solution must hold with 2 measured from the closed end of the 
tube (in order that w and ¢ shall vanish there). This solution will now be 
determined, including the value of t, for which it holds. For smaller values 
of t, the solution holds without change if x is measured from the point at 
which flow ceases (Fig. 3). 

Now the boundary conditions and equations (19) and the first of equa- 
tions (17) are satisfied by 


8 9 4— 38 ., 
ha ge FO a Fr), 


) 3) 


> (31) 
. l 6B? 
u 4B2(1—6r?+ 9r4— 4r6) — a, { 168 }eto 374+ 276) , 
v4 o 


The parameters 8 (not to be confused with the B of section 4, which will not 
occur again) and t, are then determined from the second of equations (17), 
which becomes 


' 168% = 320(126—51B— 28?) 
17 — 


(32) 
Co 13 B 
and from equation (18), which becomes 
1/48 g2 l a(t = I (i 16B7\*) 
ree ot 2 > + + 
o\35" | 420"\1" Go }* 120960\1" a | J 
24-+-7 8B? 
Pe: 488 4 , (33) 
120 o 
The solution of (32) and (33) is most easily obtained for o %, when 


there are three roots (indeed the equation in terms of B is a cubic), namely 
B = 2-091, (, =: Zak; B 4-017, t, = 5006; 
B 25°705, ty 949. (34) 
Of these only the first corresponds to a similarity solution achievable in 
practice. The evidence for this is somewhat complex, but note first that 
the other two solutions have u > 0 on the axis of the tube as well as near 
the walls, with therefore a region of negative wu somewhere in between. 
Thus the velocity profile has at least three ‘waves’ in it, as against two 
for the first profile. 
Now in reality there are probably infinitely many similarity solutions 
of equations (10) and (11) under the boundary conditions (12), with velocity 
profiles having every number of waves in them from 2 to oo. The method 


of obtaining approximate solutions which has been used is only suitable 
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for profiles with very few waves, say 2 or 3. Hence only these emerge from 
it. (The third solution already is probably hopelessly inaccurate.) 

The existence of this infinity of similarity solutions is typical of almost 
any physical problem under homogeneous boundary conditions like (12). 
(Mathematically we have an eigenvalue problem.) In problems with 
damping, the one corresponding to reality is the one with least damping, 
which here is evidently the one with fewest waves (since the damping is 
viscous). For a linear problem the general solution would be a linear 


+ and the steady solution 


combination of the different similarity solutions, 
is that similarity solution which damps out least rapidly, being attained 
when the others have become small compared with it. For the present 
problem, which is non-linear, the solutions other than the similarity solu- 
tions are related in a more complicated way (see section 6), but it is reason- 
ible to suppose that of the similarity solutions it is still only the solution 
with fewest waves which is reached as a steady state in practice. This 
conclusion is supported by the plausibility of the resulting complete (N,, t,) 
curve (Fig. 7). 

Since only the positive root of (33) is relevant, it is the only root whose 


variation with o need be considered. When this is calculated the variation 


of 8 with o is found to be quite negligible; t,, however, varies considerably. 
For o 1, for example, the root is 

B—2-099. t, = 227. (35) 
This value of t, is 27 per cent. below the value for o 0. 


Since by (15) and (31) 


¥, 12—4f 6—28 
\ | (1.2 “) kn, (36) 
™ 5 5 
the Nusselt numbers for o © anda | respectively, by (34) and (35), are 
N 0-364 (oc 0, t, 311). N,, 0-360 (oa l, t, 227). (37) 


These points can be marked in on the (N,,,t,) curve (Fig. 7). 


Note that, whereas in the present extreme case of the problem ¢, is over- 


estimated (for given NV) by putting o x1, the reverse is true in the other 
extreme case, as (22) shows. For some intermediate case, the (N,,,¢,) curves 
for o © and o 1 must cross, and in this region the approximation 
obtained by putting o © must be particularly good. 


For smaller values of t, than that (say t,) which is characteristic of the 
similarity solution when it fills the whole tube, ¢t and u are given by equations 
(31) with t, replaced by t, and « by a—lU(1—t,/t,). For 2 < U(1—t,/t,) there 


For an example of this in heat-transfer theory, consider the temperature distribution 


in laminar flow through pipes (5, p. 616). 
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is no motion and no variation of temperature (Fig. 3), and the top of this 
region plays exactly the same role relative to the similarity flow above it 
as does the closed end of the tube in the case discussed above. The Nusselt 
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Fic. 4. Velocity and temperature pro- 
files for the laminar similarity solution. 


number is therefore unchanged if in it / be replaced by the ‘effective length’ 

of tube /t,/t,. If this is not done, therefore, we have 

N.=N,,2 fort, <t, (38) 
t, j 

where the suffix s refers to the similarity values (37). On the (N,,, t,) diagram 

(Fig. 7) the relation (38) gives straight lines pointing downwards and to the 

left from the points representing the similarity solutions. 

It is interesting that the velocity distribution (measured relative to «//a’, 
see (9)) is practically unaffected by the value of o (since by (31) and (32) 
it depends only on £). Only the axial temperature gradient (proportional 
to t,) is so affected, owing to the additional thermal buoyancy effect needed 
to overcome the inertia of the fluid. The velocity and temperature profiles 
are shown in Fig. 4. The negative velocity in the middle is greater than the 
positive velocity at the outside because the larger values of r correspond to 
a larger area of tube. Similarly the temperature gradient increases at first, 
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as r decreases from 1, because it is causing heat transfer over a decreasing 
area at a rate initially changing little. 


6. Laminar flow filling the tube without similarity 

We shall now seek to satisfy the equations (17), (18) and (19) in the 
limiting case o = 00 by means of velocity and temperature profiles which 
fill the tube but whose shapes vary with x. The form of the profiles is very 
similar to (31). In fact we take 


5.2 , Sp—! 4—3f . 
t = t,.1—Pr?+ p- ae pai if 76), (39) 
+) a) 
481, . 2 ; -_ 
u (1 —6r?-+ 9r4— 4r6)— Jt. (r?— 3r4+ 27’). (40) 
dt ,/dx 7 
These satisfy the boundary conditions on r = 1 and equations (19) and the 


first of equations (17). In (39) and (40) both ¢, and B are functions of 2. 
Substituting (39) and (40) in equation (18), with o = o, we find that 


dt, 57608 


' (41 
dx 24+ 7B 
while the second of equations (17) becomes 
17 ,450+-938-+ 14B2 12—4 
ie p+ 148 Seo fl (42) 


dx 604800 5 


Eliminating x between (41) and (42), an equation relating t, to B results, 
| 5 “e 
namely 1 dt. 18(93-+ 288) 
t,dB  3024—576B—387B?— 14p3" 


The cubic in (43) has the roots given by the values of 8 in (34) for which 


(43) 


. similarity solution (with 8 = constant) exists. The solution of (43) in- 


volves these roots in the form 


t, = const(2-091 — B)—°867(4-017 + B)+0211(25-705 + B)—0 9544, (44) 
Since ¢, must be 0 at the closed end x = 0, this value of 2 corresponds to 
8 = —4-017. As a and t, increase 8 must increase towards the value 2-091 


corresponding to the similarity profiles. 
Hence, by (41), 


- dt (B), t te(B1) (45) 


57608 
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if 8, is the value of B at the orifice x = 1, and by (15) 


Bi 








12—48 244-78 
; t(B) 5 57608 dt.(B) 
N, ——  eOI7 ' : (46) 
Bi 
> 244-78 
t(B,) | arog Ue) 
4-017 


These forms of t, and V, have been written in such a form as to be functions 
of 8, independent of the constant in (44). They were computed by numerical 
integration, in which special account had to be taken of the singularities of 
the integrand. The adjacent shortened table indicates their behaviour. 
The peculiarities near the top of the table are due to the zeros of (43) and 
of the integrand of (45) at B -3-321 and B —3-429 respectively. 
A more extended table was used to plot the graph of N,, against ¢, shown 
in Fig. 7 (in which the looping of the curve corresponding to the above- 
mentioned peculiarities is invisible to the naked eye). This graph runs 
continuously into the point t, = 311, N, = 0-364 given by the similarity 
solution for ¢ = «. This limiting case occurs as B, > 2-091, when almost 
all the change in x (see (45)) occurs for 8 near 2-091, so that the similarity 
profile fills almost all the tube. 
TABLE 2 


. . 

"1 1 Na 
f'O17 5620 | 2°51 
3°207 5950 2°77 
2°517 S500 ol 
1°767 5310 59 
I'O17 4540 2°96 
0°207 4210 7 
0493 3440 O4 
1*233 2500 2°45 
I 3 1351 1°54 





But in all cases of the present solution the quite different velocity and 
temperature profiles represented by 8 4-017 are attained at the closed 
end x = 0. These profiles are illustrated in Fig. 5. The rising velocity and 
local maximum of temperature at the centre are due to convection from 
the horizontal end of the tube. The transition from this form of profile to 
one with w < 0 in the centre (which occurs for 8 > —3-429) corresponds 
to a flow like that illustrated in Fig. 6, with a standing ring vortex at the 
bottom of the tube, round which the main flow passes.t+ 


+ It is probable that something similar actually occurs in the boundary-layer flow of 
section 4, being also due to free convection from the horizontal closed end. But even if it 
does, it can have little quantitative effect on the boundary layer flow. 
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The joining up of the present solution with that of section 5 is effected 
in the next section. 
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Fic. 5. Velocity and temperature pro- Fic. 6. Sketch of ring-vor- 

files very near the closed end. tex at closed end. 
7. Synopsis of the laminar flow régimes 
The calculations of sections 4, 5, and 6 are summarized in Fig. 7. For 

t, > 3400 (if o = oc) a flow is possible in which there is a central core of 
unheated fluid. The diameter of this core is least at the orifice, but the 
volume flow of fluid is greatest there. This is because, all the way down, 
the fluid in question is being drawn into the boundary layer, where it is 
heated and forced to change direction. The logarithmic graph of N,, against 
t, in this region has the straight asymptote (corresponding to the vertical 

nd plate solution) shown, but it is noticeable that the deviation from this 

nail straight line is almost negligible even at the crossed point, where the 

snd solution ceases to have validity. This is because of the approximate 

om cancelling out of two effects. The flow of cold fluid down the centre would 

. to increase the heat transfer (forced convection) if it were not that the whole 

nds scale of the motion is reduced by the increased viscous stresses. 

the Again, for t, < 5600 (also if o = 0), a flow is possible in which the 
minimum temperature over a cross-section increases from 7; at the orifice 


to 7, at the closed end. Naturally the heat-transfer rate is substantially 
lower. (Indeed, fort, < 311, the fluid in part of the tube is at rest, or at least 
its motion is negligible.) What is perhaps even more serious from a practical 


092.24 Ee 
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point of view is its distribution along the tube; the local heat-transfer rate 
decreases steadily from the orifice to the closed end (except for a negligibiy 
slight increase just before the end associated with the ring vortex of Fig. 6), 

[t is not clear from the figure at precisely what value of t, the transition 
between the two régimes occurs, and it may be indefinite in the sense that 
there is a small range of values such that the steady state reached depends 
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on the initial conditions. The heat-transfer rate might then exhibit 
hysteresis as the temperature 7), was increased or decreased with ¢, in this 
range. But on the whole, from the form of the curves in Fig. 7, one would 
expect the transition to occur at the lower end of the range, near the 
crossed point.+ If this is so, then one might recommend the tube as a heat 
transfer device for say ¢, > 4000, and adopt the ‘fourth root’ law (22) in 
this range. 

With the aid of the limiting case curves for o = 1 one can visualize how 
the general shape of the intermediate curves discussed above might be 
altered for varying oc. Generally speaking, conditions change linearly with 
respect to o-!. 

The rough rule emerging from the foregoing, that there is a critical value, 
equal to about 4000, of t, = G,(a/l) = G(a/l)*, above which _N, 


a 


varies as 
+ It is interesting to inquire what property of the non-similarity flow characterizes this 
point at which the transition to the boundary layer régime is made. Actually the most 
striking such property is that about half the tube has a temperature profile with a local 
maximum at the axis, and half (the upper half) has one with a minimum there. Thus, the 
non-similarity solution persists until the average profile is just about to cease to have a 
minimum at the axis, and then the transition is made from this profile (already considerably 
flat at the axis) to a profile with a completely flat portion in the middle (outside the boundary 


layer). This observation will be found to be of suggestive value in the turbulent theory 
(section 11), 
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tt, 


example if all quantities are fixed, including / but not a, then there is a 


and below which it varies as ¢,, is evidently of great significance. For 


critical radius a (for laminar flow), equal to about 1(4000/G;,)', such that 
for greater radii the total heat transferred is proportional to a, and for 
smaller radii to at. Now one can imagine applications in which the number 
of holes used would depend on the radius, and might perhaps vary somewhat 
like a~? (which would keep the porosity constant). Then the critical radius 
mentioned above wouid be the radius for maximum heat transfer. Actually, 
of course, a larger value than 4,000 in the formula for critical radius would 
preferably be used, owing to the approximate nature of these calculations, 
until they had been more closely checked experimentally. 

\gain, for fluid with given p, c,, and o, the heat transfer will increase 
with viscosity in the upper range of values of ¢t, (because t, oc v-* and so 
Q oc kN, pc, a 'vN, x vt} oc v') and decrease with viscosity (like v~!) in 
the lower range. This means that there is an optimum viscosity, given 
everything else. It also means (what is more important) that the general 
effect of turbulence (which to the roughest approximation is like an increase 
of v without change of p, c,, or o) will be to decrease heat transfer in the 
lower range of t, from what it would be in laminar flow, and increase it in 
the upper range. Physically, this is because in the lower range the impeding 
of the flow is more important than the increased efficacy of the flow in 
transferring heat, and in the upper range the relative importance of these 
effects of turbulence is reversed. The quantitative extent of this distortion 
of the curve will be seen in section 11. 


8. Turbulent flow: general considerations 

There are two factors, in the laminar flow which has been described, 
which tend to make it become turbulent. First, the cross-sectional distribu- 
tion of shear has a maximum in the fluid. Secondly, the temperature 
gradient upwards is negative. Also there are two factors resisting transition 
to turbulence, namely, viscous damping and the proximity of the walls. 
But these latter factors are present in straight pipe flow, where there is no 
shear maximum (except at the walls) and no temperature gradient. Hence 
in the present problem the transition Reynolds number may be expected 
to be considerably less than in straight pipe flow, and indeed Hermann (9) 
found that for free convection from both vertical plates and horizontal 
circular cylinders transition occurred when a Reynolds number based on 
maximum velocity and on boundary-layer thickness was about 300. 

Wall roughness is unlikely to be important, because the velocities are 
only moderate in free convection flow. According to the predictions 


below, even for the largest modified Grashof numbers contemplated in the 
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engineering application, the velocities are such that walls on which there 
were no pimples of height exceeding | per cent. of the pipe diameter would 
be aerodynamically smooth. But conditions at entry to the tube will be 
considerably disturbed, whether the fluid enters from a large closed reservoir 
(continually stirred up by the discharge of heated fluid), or from a circulating 
system, or (as in the system of Fig. 1) from a region of mixing of the hot and 
cold fluids which have come together at the cireumference and then flow 
away again down the middle. Thus it is advisable to assume low values of 
transition Reynolds numbers, corresponding to large entry disturbances, 
The existence of three régimes of flow, corresponding to those of section 2 
in laminar flow, and appearing in different ranges of values of G, and I/a, 
is to be expected. But the range of //a corresponding to the first régime, 
with the boundary layer not filling the tube, is doubtless reduced con- 
siderably owing to the more rapid thickening of turbulent layers. Condi- 
tions in the limiting case of this régime, when it is equivalent to convection 
from a vertical plate, can be obtained from Saunders’s experiments (2) for 
G, < 10", and by extrapolation for values of G, beyond this. These yield 


the empirical formulat - » 
N, = 0-1143. (47) 


( 


3ut also, by analogy with the results of section 4, one may suppose that 
this is a good approximation throughout the limited range in which the 
first régime obtains, since the additional forcing of the convection by flow 
down the centre is balanced by a general reduction of velocity due to 
increased relative shear. But the cessation of this régime when the flux 
outside the boundary layer ceases to be a maximum at the orifice is made 
even more certain now that any fluid which has ever been in the boundary 
layer must continue to be turbulent. 

It is the second and third régimes (which, as is clear from the above, have 
additional importance in turbulent flow) that present the most novel 
features to the fluid motion scientist and call for the formulation of new 
hypotheses if they are to be elucidated. The turbulence fills the whole tube. 
but the velocity distribution is very different from that of turbulent pipe 
flow, since it changes sign in the middle. 

In these régimes, doubtless, the general scale of both the mean motion 
and the turbulence will decrease between the orifice and the closed end. 
just as it does for laminar flow. Probably, indeed, there is a short length 

+ Saunders gives this formula with both N and G@ based on the plate height /, but since 
] appears to the same (first) power on both sides of the equation it may be replaced by 
any other length, here the radius a. There is, of course, no a pricri reason why heat transfer 


per unit area shall not be independent of the plate height. For, higher up the plate, the 
heat-transfer rate is increased owing to the higher velocities reached and decreased owing 


to the larger boundary-layer thickness, and these may cancel out 
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of tube near the closed end where the local Reynolds number is so small 
(say about 100) that turbulence is rapidly damped out and the flow is 
effectively laminar. At least it will be convenient, in the theory to be 
developed, to take this transition to laminar flow as a sort of boundary 
condition at the closed end. 

Now before developing any special hypotheses regarding turbulent 
exchange coefficients, we can at once observe an important difference 
between the properties of the similarity solution in laminar and turbulent 
flow, by writing down the equations of the mean motion and remembering 
that broadly the turbulent shearing stress varies as the square of a typical 
velocity. We may write the equations of the mean motion (compare (3) 
and (4)) in the form 


cu oV V i a oT l ¢ 
0, l 1 V Rq), 48 
aX 'oR'R >X*' dR RoR” ”? _ 
oU oU . | eo 1 ¢ cp 
f ; J ] P | Rr), 0. 49 
( X ¢ R ; Pp rf X R C R' Cc R 


Here l’, V, 7, p signify local mean values of the quantities so denoted 
previously; q signifies the mean rate of radial heat flow per unit area per 
unit time, divided by pc., (the heat capacity per unit volume); and 7 signifies 
the mean shearing stress divided by the density p. 

To eliminate the unknown pressure gradient @p/¢eX in (49), the device 
used in section 3 (subtracting the equation in the form it takes at the wall) 
is not really suitable because the gradient of shear stress at the wall is not 
i significant quantity in turbulent flow. It is better to subtract the mean 
value of the equation over the whole cross-sectional area (obtained by 
multiplying by R, integrating from R = 0 to R = a, and dividing by $a’), 


giving 
.oU ~0U - a ; : . a 2r 
l J RU2dR (T—T.)- Rr)+-—®, (50 
c X c R a” Cc zx R ; u m) R Cc R| ; ) a \ ; ) 


0 
where 7, is the skin friction (namely minus the value of 7 at the wall) and 
T,, is the mean temperature across a cross-section, namely 
2 . 
7 _ | RT AR. (51) 
a* 


0 
It might as well be said at once that the left-hand side of (50), repre- 
senting the inertia of the mean flow, will always be neglected compared 
with the right-hand side, which represents the effect of buoyancy and 
shearing stress. This was found to be a reasonably good approximation in 
laminar flow, and is presumably better still in turbulent flow when shearing 
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stress is relatively more important. Thus the final form of the momentum 
equation, which with (48) governs the motion, is 

a. 
1 < 279 


(Rr) + 


i. — (), 52 
“ROR ; — 


Now in turbulent flow, to a very rough approximation, the order of 
magnitude of 7 at a cross-section varies as the square of the order of magni- 
tude of U, while that of g varies as the product of those of U and 7. It 
follows, by (52), that in a similarity solution the scale of 7’—7) must vary 
with X like the square of the scale of U’; and by (48) that the scale of 67'/eX 
must be the same as that of 7’. This means that the variation with X is 
exponential, not linear as with laminar flows. In fact to the rough approxi- 
mation mentioned 7'— 7, varies like eX and U like e!4*, for some A. 

More exactly, there are two departures from the truth of this conclusion. 
First, since the simple laws for shearing stress and heat flow are only 
approximate, A varies somewhat down the tube. In other words, if we 
replace €7T/0X and eU/eX in the equations by A(7—T),) and 4AU, and 
solve them (see section 10), we shall find that A varies. Secondly, the 
exponential mode of variation must cease somewhere or the tube would be 
infinitely long, and evidently this occurs at the transition to the laminar 
region, in which the variation with x becomes linear. 

After a hypothesis regarding the turbulent exchange has been formu- 
lated in section 9, the similarity solution will be calculated in section 10. 
This will be found to hold for a particular value of I/a, which is a function 
of G, but not as in laminar flow directly proportional to it. For larger 
values of l/a the additional length is expected to be stagnant. For smaller 
values of //a a solution without similarity properties will be obtained in 
section 11. 

As in the case of laminar flow the equations will be used only (i) at 
particular points and (ii) in an integrated form. The integrated forms of 
equations (48) are 


a a 
oo e fe ra i 
RUdR 0, . RUT dR = aq, (53) 
aX . 

0 0 
where q, is the value of ¢ at the wall. On the other hand, this integration 
gives no additional information in the case of equation (52), as is clear from 
the manner in which it was obtained. The equations of heat and momen- 


tum on the axis take the form 


»~ oF eq 9 Or 27 i 
us 2 | a he —2=6. (8 
| Ee 0 ( 2), 0 - m)R-0 ( ?) a . 
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They will also be applied at the outer boundary of the so-called laminar and 


transition sub-layers (see section 9). 


9, Hypothesis concerning exchange coefficients in confined turbu- 
lence 
It is necessary to have a theory of the influence of the turbulence on the 
mean flow sufficiently detailed for skin friction and heat-transfer rate to be 
predicted; the kind of theory which assumes the values of these and then uses 
the integrated equations of motion to determine other properties of flow is 
useless here because of the depth of our ignorance. Now detailed theories 
of the influence of turbulence on the mean flow have always been based on 
hypotheses concerning the distribution of the ‘exchange coefficient’ or 
‘eddy kinematic v iscosity’ e defined by the equation 
eu eu 
“@R'' OR’ 


Here the turbulent component ¢ ¢U’/¢R is less than the laminar component 


(55) 


veU eR only near the wall; but there it becomes very small indeed. 
This is clear from Reynolds's expression for it as u’v’ (the covariance of 
ixial and radial flow velocities) which, as Reichardt (6) points out, falls 
to zero like (a— R)%, since u’ is proportional to a—R and v’ to (a— R)?. 
\way from the wall e greatly exceeds v. In this region it was for many 
years usual to apply mixing-length theories, which relate e to the local mean 
shear. These are now discredited and their success in certain problems 
scribed to particular coincidental circumstances. In the present problem, 
where the velocity profile has three stationary values, the exchange coeffi- 
cient would on momentum transfer theory have three corresponding zeros, 
which is absurd. Far from the turbulence being so small in scale that it 
reaches an equilibrium intensity depending only on local flow parameters, 
it rather seeks to fill the whole tube, at the walls of which the eddies have 
been described (7) as being ‘anchored’. Near the wall the size and intensity 
of those eddies responsible for exchange is limited by the proximity of the 
wall, and it appears likely that this rather than the local mean shear is 
responsible for the well-known logarithmic profile. 

These considerations lead to the hypothesis that the distribution of 
exchange coefficient in confined turbulent flows may depend solely on 
position and on some factor representing the scale of the turbulence; a 
reasonable such factor is the shear stress it produces at the wall. According 
to this hypothesis the distribution across a cross-section of factors tending 
to create turbulence is of no importance; only the average is of importance, 
as determining the scale of turbulent exchange. This is reasonable if turbu- 
lent energy is created in the form of large eddies. The mechanism by which 
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the turbulence transfers quantities in proportion both to their local gradient, 
and to a function characteristic of the local scale of the smaller eddies (the 
exchange coefficient), could then be that suggested by Townsend (8): the 
smaller eddies are rendered more nearly isotropic by energy transfer down 
the scale of eddy sizes, and at the same time oriented by the local shear 
(so that they cause a shearing stress), and a balance between these two 
effects determines the exchange coefficient. 

Now in one problem of confined turbulence with smooth walls (flow 
through a straight pipe under a pressure gradient) the exchange coefficient 
has been determined from measurements of velocity distribution (with the 
help of theoretical considerations), for a wide range of values of the skin 
friction, by Reichardt (6). This can be taken over, according to the above 
hypothesis, for any flow in a circular pipe, including that at present 
discussed. 


Reichardt’s argument will be briefly recapitulated. In the flow through 


Ly) _R 
dR be 


Now it is known that outside the ‘laminar’ and ‘transition’ sub-layers the 


a straight pipe du 


(e (56) 


velocity distribution is given to a good approximation by the formula 


; \T a—R , = 
l log (Avro ) A204, Ax. (57) 
K v 
Reichardt finds departure from this in the experimental values only very 
near R = 0, where the velocity gradient does not change sign quite dis- 
continuously. Now, in the region outside the sub-layers, « > v, so that 
equations (56) and (57) give 
: Ria—R 
€ K NT) ) (58) 
a 
In particular when (a— R)/a is small (near the wall), we have approxi- 
mately , . 
. € kK NT (a R), (59) 
corresponding to the profile (57) holding in the region (sufficiently near the 
wall but outside the sub-layers) in which the shearing stress is effectively 
constant, i.e. in which , 
dl 


hie?) am 


This profile has been verified for a variety of flows near a wall. 


To: (59a) 


However, at the wall itself « must fall off like (a— R)’, as observed above. 


And, indeed, if equation (59) held right up to the wall, the solution of 
equation (59a) with U 0 at R = a would be of the form (57) outside the 
sub-layers but with A K instead of A > 20K. Reichardt finds that a 
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ient. simple smooth curve for ¢, which starts to grow from the wall very slowly 
(the like (a— R)® and then becomes asymptotic to (59) at about 
- the 


(a— R)v7,/v = 20 
own 


hear (corresponding to the boundary of the sub-layers), gives a formula (57) for 


two | U outside the sub-layers with the correct value of A. This gives a simple 


view of the laminar sub-layer as the region in which « is negligible and the 
flow transition sub-layer as the region in which it is changing from propor- 
‘ient tionality to (a— R)* to proportionality to (a— R). 
| the Now it is important to notice that when this theory is applied the exact 
skin | form of ein the sub-layers need not be used. For, provided that the shearing 
ove stress is effectively constant in the sub-layers, the velocity distribution in 
sent them must be the solution of (56) which vanishes at the wall. In particular, 

just outside the sub-layers it must take the form (57). The proviso is 
ugh reasonable because even if some lateral gradient of shearing stress is 





required in the sub-layers to balance other axial forces, the sub-layers are 





(56) so thin that the shearing stress in them cannot vary much from 7». 
Hence a satisfactory treatment of problems of flow in pipes should result 
the from using a form of « appropriate to flow outside the sub-layers (and 
;, neglecting viscous stress there), and applying the wall-boundary condition 
7 ( ‘70 log| Avr9— =| ~O as Roa. (60) 
kK vy 
very 
die- The velocity distribution so determined will depart from the true one only 
that in the sub-layers, where its detailed distribution has no effect on integrated 
quantities such as those in (53). The boundary condition (60) corresponds, 
therefore, to expressing the correctness of the momentum equation in the 
(58 sub-layers, instead of at the wall as in the laminar theory. 
[t will be assumed that the same exchange coefficient governs the flow 
_ of heat as of axial momentum (the Reynolds analogy). This reasonable 
(59) hypothesis has never been disproved. Note that it does not (as is sometimes 
the believed) imply the equality of skin friction and heat transfer from the wall 
vale (when measured in suitable units), even for o = 1, but rather causes depar- 
, tures from equality in the right direction (as Taylor (10) showed). On the 
9 a) hypothesis of 


q = (e+) (61) 


eR 


solving the equation gq = qg, in the sub-layers. The solution will take a 


ve. | The boundary condition on 7’ corresponding to (60) can be determined by 
| particularly simple form, 7'— 7) (qo/T))U, ifv = x, i.e. ifo = 1. Since 


dependence on o, for a given value of the modified Grashof number, is not 
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expected to be very great (judging from the results for laminar flow) the 


assumption o = 1 will be made henceforth. The boundary condition is then 
mem, ! a—k , 
T—T,- fo log( Avro >~O asR-a. (62) 
Kv, v 


It is worth remarking that for any o (62) would still be true if the value of A 
were modified; and that the true A lies between A and Aco.t 

When the boundary conditions (60) and (62) are used the only remaining 
problem is the choice of ¢ in the region outside the sub-layers. As remarked 
above, Reichardt finds that the equation (58), which with the boundary 
condition (60) is equivalent to the velocity distribution (57), produces 
velocities differing slightly from those observed; in fact the discontinuity 
in velocity gradient is slightly smoothed off in the observed distribution. 
This corresponds to a value of « which does not fall to zero at R = 0 but 
rather flattens out to a central value of about half the maximum of (58), 

However, in the present problem we are not interested in the fine struc- 
ture of the velocity distribution. In fact it is necessary that we neglect it, 
since to apply a Pohlhausen-type method of calculation it is necessary to 
use velocity profiles of simple analytic character, whose shape does not 
suddenly alter near the centre. Also, since the centre is a point at which 
the equations will be forced to hold, it is desirable that it should be made 
typical of neighbouring points. 

Hence, finally, the equations of section 8 will be used with 


R(a— R) eu : R(a— R) eT 
q Kn 


: om (63) 
a cR a oR 


t= Kv, 
and the boundary conditions (60) and (62) at the wall. Here, +, and q, are 
functions of x, which are ultimately to be determined. Because the exchange 
coefficient vanishes at the centre the condition of evenness of U’, 7’ (in the 
sense that 6U/éR and éT/éR vanish at R — 0) has to be abandoned, being 
replaced by a condition that U and 7 are finite at R = 0, which under this 
circumstance is just as restrictive as is the condition of evenness for non-zero 
exchange coefficient. 
The equations (53) and (54), with (63), can be put into non-dimensional 
form by the substitutions 


, v a y* . a 
l -U, T—T, = —-—.t, X == ¥, R = ar, (64) 
ak vfas K2 
+ The Grashof numbers and modified Grashof numbers coincide for o 1. Hence the 


results will be given in terms of Grashof number. But it remains true that for fluids with 
o 1 the results probably have slightly more value when the Grashof number is interpreted 


as a modified Grashof number, and the above remark provides a reason for this. 














whe! 


with 


hold 


loca 


Wit 


An 
em] 


10. 

A 
is 0 
(66 
whi 
slov 
apy 
mu 

A 


whi 


ant 


tut 


ant 








the 


| 
nen 


ary 
ices 


ty 


on. 














FREE CONVECTION IN TUBES 


when they become 


1 1 
Cc > 
ru dr = 0, | rut dr = —QR, (65) 
mn Cx Fp 
0 0 
1 
: ct ot ¥ cu P _— 
with u a, tj | rt dr = 2R—-+ 2R? (66) 
Cx Cc} ° d cr 
0 
holding on r = 0. Here R is av7,/v, the Reynolds number based on the 
local wall shearing stress, and 
a} 
v/arc = 
Q J 10 : (67) 


Kv?v7, 
With the new variables the boundary conditions (60) and (62) become 
u—Rlog{ AR(1—r)} > 0, t—QlogfAR(1l—r)} > 0 asr—>1. (68) 


An advantage of the equations in the form (65), (66), and (68) is that the 


empirical constant K has been eliminated. 


10. The similarity régime (turbulent flow) 
As stated in section 8, the similarity régime with T’—T, « e*, U ac e!4* 
is only approximately achieved in turbulent flow. However, if in (65) and 


by 2Aut and ét/éx by At, then we find solutions for 


66) we replace c(ut)/éa 
which both the shape of the velocity and temperature profiles vary only 
slowly with x. The substitutions are therefore permissible, similarity being 
approximately achieved in the sense that the profiles change shape with x 
much more slowly than they change scale. 

\s profiles satisfying (68) we choose 
u = R{log(l—r)+B—C(1—r)], ¢t Qi log(1—r)+ B+ D(1—r)], (69) 
where B log(AR). Then the first of equations (65) gives 
C = 3B-3, (70) 
ind the second of equations (66) gives then 
Q(2+2D) = 2R2(3B—%8). (71) 
The first of equations (66), with d/dx replaced by A, gives (after substi- 


tuting for C from (70)) 


R(3—2B)AQ( B+ D) 2RQ(1+- D) (72) 
and the second of equations (65), with d/dx replaced by 3), gives 
SARQ(2—4B+1D—}BD) OR. (73) 


Eliminating A from (72) and (73) we obtain 


(81—24B+-17D—6BD)(D+1) 12;(9—4B)(B+D). (74) 
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The solution of (74) for D is given as a function of B in Table 3. It 
is interesting that in the range of values of B which is of greatest 
practical interest (which will be found to be the range 5 < B < 8) the 


coefficient D is practically constant at 10-5. The associated values of 


R = e”/A (with A taken as 9), Q from (71), and t, = Q( B+ D), the value 
of t at r = 0, which is a sort of local Grashof number, are also given. It is 
worth noticing that t, is approximately 10°-1 in each case. The product QR 
is a non-dimensional form 

vga 


K ys Yo 


of the local heat transfer rate. Ny is the Nusselt number based on the total 


OR 


(75) 


heat transfer from the part of the tube below the cross-section X, and on 
the length X. Thus 





xX x 
v (pt, de 2a dX Daa! ie 
iV 7" yr | PCy 0 HE Oe a) Yo G4 
i(T,—T)2aXa/Xx | ° (T_T) | 
0 0 
x 
UF ORd QR(; BD—13D+13B—2) (78) 
— dx va = : ; 76 
Kt, , KQ(B+D) 
0 
[t has been computed below taking AK = 0-4. 
TABLE 3 
B\ D R Q t. QR Ny A 
4] 12°5 6:07 56 | 99.108 | 34.102 | 131 | -468 
5 | 106 | 16: 66.10? | 1:04.10 | 11.104 | 7:8 271 
6 | 10°3 | 44°8  ! G5.10* | 1:05.10° | 2°9.10° | 24-7 “185 
7 | 10°4 122 | 58.104 | rror.10° | 71.10% | 84°3 138 
8 | 10°7 331 5°0.10° | 9°3.10° | 1°64.10% | 272 | -109 
g | 111 goo 4°1.10° | 82.107} 3:7.10°| 8 


550 "059 








Finally A (obtained from (72) or (73)) tells us how the quantities vary with 2, 
in virtue of the relations 
dt, * d(log t,) 


kis x ‘ 77 
- A, 2 X (77) 


The integral (77) for 2 is made definite by considering the transition to 
laminar flow, which, as explained in section 8, must occur near the closed 
end when the Reynolds number falls to such a value that disturbances are 
rapidly damped out. If the transition value of Bis B,,, then by a combined 


use of (77) and the result of section 5 for o 1 we have 
B 
- Bw * d(log 
2 K2lel r) | (logt.) (78) 
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bearing in mind that ¢,(B,) is the Grashof number based on the radius and 
on the central temperature at the transition point, and that x = k?X/a. 

The value of By, is rather difficult to forecast at all precisely, but fortu- 
nately equation (78) is very little sensitive to it. It is almost certainly 
between 4 and 5. Thus for B = 4 the Reynolds numbers based on the 
radius and on the maximum upward and downward velocities respec- 
tively are 15 and 50; while for B = 5 they are 70 and 250. Taking into 
account Hermann’s observation of 300 as a typical Reynolds number 
based on maximum upward velocity and on boundary-layer thickness at 
which turbulence set in for free convection flow up a vertical plate, but 
remembering that a much lower value is required if existing large distur- 
bances are to be speedily damped out, and remembering also that the net 
difference of velocity across the inner shear layer gives a Reynolds number 
obtained by adding those mentioned—but, on the other hand, that the 
complete confinement of the flow by walls must assist the damping-out 
process—we may reasonably fix ona value B, = 4-5, at which (by logarith- 
mic interpolation in the above table) 

R= 100, Q= 200, ¢,= 32.16, QR = 2000, N, = 3-46, 

A = 0°340. (79) 
Then the length of tube with laminar flow is 3-2.10%/227 = 14 radii. This 
assumes that the similarity régime occurs also in the laminar region. This 
can be checked by comparing the value of N, in the laminar and turbulent 
flows at the join. This is the only parameter (representing total rate of heat 
transfer across the join) which must remain continuous however compli- 
cated the mixing process at the transition point. But for the laminar 
similarity solution with //a 14 and o = 1 we have 

N, = 14N, 14(0-360) = 5-0, 

which is of the right order of magnitude. (Actually any of the non-similarity 
flows filling the tube give also the right order of magnitude for N,, and 
some of them give much smaller lengths of tube with laminar flow; it is 
quite likely that one of these actually occurs. But in any case, even the 
14 radii are small compared with the total length of tube required for the 
similarity régime.) 

It should be noted that the whole theoretical basis of the method of 
calculation is at fault for values of R as low as 10. In fact for B < 6 
(R < 45) the assumption that the sub-layers fill only a small part of the 
tube is incorrect; and at the transition value B = 4-5 the laminar sub-layer 
(a—R)v7)/v < 10 actually fills the whole tube. Thus for B < 6 the 
formulae can be regarded only as interpolations between the laminar and 


fully developed turbulent flow régimes. 
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The integration in (78) was carried out graphically, giving a relationship 
between X/a and B. This is also the relationship between //a and the 
orifice value of B, which determines, by Table 3, the orifice values of t, and 
N,, which are G, and \,. The results are finally obtained in the form of 
and G, = (l/a)?G, and N, and N, = N,a/l against //a for the 
turbulent similarity régime. The results were also determined by an 


graphs of G, 
approximate method which will be developed in section 11, and there was 
good agreement. These results are shown in Fig. 9 below as the curves 
marked D 10:5. It is seen that for the larger Grashof numbers the length 
of tube required to attain the similarity solution is very great, and also that 
the heat transfer remains only very moderate, as in the laminar similarity 
regime. 


11. Turbulent flow filling the tube without similarity 
Now also in the non-similarity régime we shall seek solutions of (65) and 
(66) in the form (69), and then the condition of conservation of mass and the 
momentum equation at the axis again give (70) and (71). But in the other 
two conditions d/dx can no longer be replaced by a constant, and accordingly 
they are 
(4B 9) © (Q(B D)| 4Q(1+ D), (80) 
and 


5 [RO(6B 17)D+-3(8B—27)}| = 72RQ, (81) 


instead of (72) and (73). 

Now it would be quite possible to solve (80) and (81) numerically, subject 
to B = log(AR) and (71), as simultaneous differential equations. However, 
it would be very laborious, and the author accordingly took a simpler course, 
suggested by a property of the similarity solution. This was the property 
that D was practically independent of 2. This makes it likely that, to a 
sufficient approximation, in the solutions appropriate for tubes which are 
shorter (for a given Grashof number), D is also a constant, but not the 
same one. This assumption makes it possible to discard one of the equations 
(80) and (81), and it is reasonable to discard the former, which expresses 
conservation of heat at the axis (and is in any case somewhat doubtful 


because of uncertainty regarding the form of ¢ on the axis), in favour of 


the latter, which expresses a net balance of heat transfer for a whole cross- 
section. Then three equations have been retained, expressing conservation 
of heat, mass, and momentum respectively. The additional parameter D 
(a constant, but different for different solutions) permits the application 
of the theory to tubes of various length—radius ratios. 
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Now by (71), on this hypothesis, 


1 dQ 2 dR 2 dB 4(B—1)dB 


(82) 
Qdx Rdx'2B—3dx  2B—3 dx 
since log(AR) B. Hence, by (81), 
dB i(—%) 6D-+24 
dx\ 2B-—3 6B—17)D+-3(8B—27) 
a2 
= C ee =p (83) 
(6B—17)D+3(8B—27) 
and so =e x,(B)- Dx'(B), (84) 
ap 
where one finds after some reduction that 
v,(B) l Be SB > log(2B—3)- const | (85) 


r,(B) = 1B2—13B—! log(2B—3)-+const} 


The constants are chosen to make a, and 2, vanish for the ‘transition’ 


value B= B, 1-5. The functions so obtained are tabulated in Table 4. 


TABLE 4 


This table shows clearly the influence of varying D. The similarity 
solution with large D = 10-5 (corresponding, by (69), to a highly peaked 
temperature distribution) is only possible for a very long tube, owing to the 
term Dx,(B) in x As D decreases we have solutions corresponding to 
shorter and shorter tubes, with the temperature distribution flatter and 
flatter in the middle. This is illustrated in Fig. 8. The limit of this state 
of affairs is given by D |, when the temperature gradient at the middle 
vanishes, and therefore a transition to a flow régime with a boundary layer 
not filling the tube may reasonably be expected. (Note that for D - l 
the axis would become a local maximum of temperature.) This is precisely 
analogous to the state of affairs exhibited in the laminar case (section 6), 
when the transition from the non-similarity flow filling the tube to the 
boundary layer type of flow occurred (section 7) very near the value at 
which the average temperature distribution ceased to have a minimum 


value on the axis. 
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For each of the possible flow régimes in the range —1 < D < 10-5 we 
need to know not only 2 but also 


1 8e?4(2 B—3)(B+D) 
























































t, = Q(B+ D) = — (86 
) A(9+4D) 
, e3( BD—4D+i1B-2 
and Nx (is 1 = s) (87) 
; Ak(B+D) 
- 
| 
12 9 
tJ 
| 
= 
s = 
WIS 
0 0:5 1 15 2 Zo 3 
logio(X/a) 
Fic. 8. The shapes Fic. 9. Variation of ¢, (Grashof number based on 
of temperature pro- radius and on temperature difference between wall 
file in the turbulent and axis) and Ny (Nusselt number based on length 


flow. of tube) with X/a, keeping (i) B equal to 5, 6, 7, 
8, 9 (ii) D equal to —1, 1, 3, 7, 10-5. 


as functions of B. From calculations of these t, and N, have been plotted 
against X/a (calculated as k-*x = 6-25x) for B = 5, 6, 7, 8, 9 and for 


D 1, 1, 3, 7, 10-5 in Fig. 9. Here no allowance has been made for the 
length of the laminar region. This is because the values of N, and t, at the 
transition value of B (namely 4-5) indicated a solution with large ¢, as most 
likely in the laminar range, so that the length at,/t, of the laminar range 
is probably negligible. In fact Ny/t, takes values which range from 9-6. 10-4 
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(for the similarity solution D = 10-5) to 2-2.10-4 (for D = —1). Regarding 
these as end values for a laminar flow, Ny/t, is the same as N,/t,, which 
takes the above values for values of t, of the order of 104, leading (since 
t. is if anything smaller than 10*) to lengths of laminar flow of less than 
one radius. 

The variation of heat transfer along the tube forag: enflowD = constant 
can be observed from Fig. 9. For the local heat-transfer rate is proportional 
to d(Nyt,)/dx. It is seen that even for D = —1 it varies very rapidly, like 
the fourth or fifth power of X/a. In a practical application the resulting 
reduced heat transfer rate near the closed end would be most unfavourable. 


12. Synopsis of the turbulent and laminar flow régimes 


In Fig. 10 a cross-plot is made from Fig. 8 for six values of l/a, after X 


has been put equal to / so that Ny = N, and t, = G,. The cross-plot is of 
NV, = N(a/l) against G,, for//a = 5, 10, 20, 50, 100, and 200, and is indicated 


by the dash-dotted lines in Fig. 10. In each case the curve is joined smoothly 
by a broken line to the laminar solution (shown as a full line), at the value 
of G, corresponding to B = 5. (A slightly higher value of B is used for the 
spontaneous appearance of turbulence, even under large entry disturbances, 
than was used above for the disappearance of fully developed turbulence. ) 
Also, at the upper end of the curve, corresponding to D = —1, a dotted 
line is shown to indicate the transition (see above) to the régime in which 
the boundary layer does not fill the whole tube. In this régime N, and G, 
are related according to equation (47). 

The position of the similarity régime is indicated on each curve by a cross. 
For //a = 5, 10 and 20 this occurs in the laminar part, but for //a = 50, 
100 and 200 there is a turbulent similarity régime instead. The dash-dotted 
curves below the crosses were obtained by supposing that the tube is filled 
by the similarity flow plus a stagnant region. None of these curves has been 
allowed to fall quite as low as the value of N,, predicted by theory for B = 5, 
for which value (as explained in section 10) the theory has little foundation; 
they have been permitted to begin departing from the theoretical curve 
from about B = 5-5, although only rising sharply at the value of G, corre- 
sponding to B = 5. 

It is seen that transition is predicted when G, is between 104 and 2.104. 
3ut it should be emphasized that this prediction, which makes the laminar 
régime play only a small part, is based on the hypothesis of large entry 
disturbances (of scale comparable with the tube diameter and typical 
velocities in the flow). With care taken to obtain smooth entry conditions 
it should be possible in the laboratory to reproduce most of Fig. 7. But in 


engineering applications the author expects Fig. 10 to be nearer the truth. 
5092.24 _£ 
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The (at first sight surprising) ‘bucket’ in each curve is due to the fact A 
explained at the end of section 7, that transition to turbulence in the régime 13, 
in which heated fluid fills the whole tube causes reduced heat transfer, 
because the impeding of the flow is then more serious than the increased 
effectiveness of the flow in transferring heat, while on the other hand the 
position is reversed once a boundary-layer régime has been set up. This | oft 
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Fic. 10. Relationship between N, 


, and G, for six values of the length 
radius ratio, with position of transition to turbulence chosen on the 
assumption of large entry disturbances. Similarity solutions shown with 


a Cross. 


the boundary-layer régime appears. The author believes that this rise has 
at any rate not been overestimated owing to the extrapolation of Saunders’s 





results far beyond the range in which they were observed, because according 





to general experience with boundary layers one would expect the variation 
of N, with G, to be as the 2th rather than the ird power, for large (G,, so that 
the upper curve in Fig. 10 may be placed too low. 

In Fig. 11 the critical Grashof number, based on length, at which the large 





rise in heat transfer occurs, is plotted against the length—radius ratio as a 
full-line curve, the broken lines giving the results of the laminar and | It 


turbulent theories separately. The transition between the curves has been tra 
made to take place between B 5 and B = 6. This curve (or, later, one the 
corrected in the light of experiment) would be used most often to determine dif 


the tube radii in terms of their length and the other constants of the 
problem. It is clear from all that has been said that to get best results the {| pa 
point (//a, G,) should lie on the left-hand side of the curve. 
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Another mode of presentation of the results, particularly useful in section 


13, is given in Fig. 12. Here a special non-dimensional form 


P92 

a xfa* 8 
N,G res Q (88) 
aTThVK 
of the heat transfer rate Q, not bringing in the tube length or temperature 
difference, is plotted against G, for various values of !/a. This displays the 
effect of length (for given radius) and temperature difference separately. 


20 on 








Logig G, 7 it) 








0 0-5 1 15 2 25 
Logio(l a) 
Fic. ll. The critical Grashof number 
(based on length), at which flow with 


boundary layer not filling tube is pre- 
dicted to set in, as a function of the 
length—radius ratio. 


It is seen that conditions in which an increase in length reduces heat 
transfer (by breaking up the flow in which the boundary layer does not fill 
the tube) are common, but that cases in which an increase in temperature 
difference has this effect are negligibly rare. 

To conclude this section it must be repeated that all the results of this 
paper, and especially those relating to turbulent flow (as being based on 
more doubtful hypotheses) must be regarded as purely tentative in the 
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absence of experimental confirmation. They probably have little more than 
suggestive value, in guiding experimental workers in what to look for. 
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Fic. 12. Curves of the non-dimensional 

form of the heat transfer rate which 

does not involve the tube length or 

the temperature difference, namely 

NG, (a fa*/27k«v)Q, against G, for 
various values of l/a. 


13. Flow in a completely closed tube, part heated and part cooled 
An attempt is now made to apply the above theory to the conditions 
illustrated in Fig. 1. Here the difficulty arises in the treatment of the 





complicated mixing process which occurs at the join between the heated 
and cooled portions of the tube, where the two circumferential flows come | 
into direct collision. It is natural to try to treat this region by supposing 
that it is of limited extent and approximating it by a discontinuous join 
between two solutions of the type discussed in sections 2-12. But if the 
flow is to be completely determined it is necessary that two conditions can 
be applied at the join, relating the flows which mix there. 

Two conditions which are evidently applicable are that the total flow 
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across a cross-section of both mass and heat are continuous at the join. But 
only the second of these is a condition that gives any information, equating 
the total rates of heat transfer in the two portions of the tube. For of 
course in all solutions the total mass flow across any cross-section is zero. 

It would be natural, as a second condition, to express continuity of total 
momentum flux at the join. But this condition need not be satisfied, 
because there is nothing to prevent a substantial pressure difference 
between the two portions of the tube from building up; and in any case the 
condition would be unsuitable since inertial forces have been neglected. 

Instead we apply a condition that the mixing is ‘complete’, in the sense 
that at least the fluid that flows away from the mixing region down the axis 
has the same temperature on both sides. This condition is probably some- 
where near the truth since the mixing is highly forced. 

The problem then reduces to determining two flows of the type of 
sections 2-12, with the same values of 7; and Q, but with different values 
in general) of all other constants. These different sets of values will be 
denoted by the symbols H for the heated portion (with walls at temperature 
T);,) and C for the cooled portion (7),). The physical constants of the fluid 
will be given values appropriate to temperatures near these values in the 
two portions. The external acceleration f (if centrifugal) may similarly 
have different average values f,;, and f, therein. Then 


Y x fh a’ 4 Al 4 Al ’ Xe ya? rn 
Gan Het ( lon T;), Goc : uf (7; -_ Toc) - (89) 


Vi Ky Vo Ko 


Hence, introducing the two known constants 


j 3 r a> 
” Ye JH @ pm mm 1H ro fo@ in 
Gan (Thu —Toc); Ta —— (Th —Toe); (90) 
Vy Ky Voke 
G Goce 
we have oF Wt. 08 ] (91) 
Gt G4, 
aH at 


as a relation between the unknowns G,,, and G,~. Expressing also that 
the heat-transfer rate 


a a. I. \ NT — n y ¢ 
Q = 2naky (Thy —T) Ny = 27ake(T,—The)Nec (92) 
is continuous, and using (89), we have 


Nic Gac rofo keveke 


“a (93) 
4 H 7 


1H ‘udu ky Yn Kn 
The product , G,, whose values in the two portions of the tube, by 
(93), have a ratio predictable in advance, is that plotted against G,, for 
different values of //a in Fig. 12. 


Now, given the constants of the problem, assume G_,,, as some suitable 





438 M. J. LIGHTHILL 


yalue about half as much as G¢;, and deduce N,, G,,, from Fig. 12 (using 
lia = l,/a). Deduce Ne Gi by (93), and infer G,¢ from the graph (using 
l/a = I,/a). Now recalculate G,,, from (91). Using this second approxima- 
tion to G,,, repeat the process. The usual methods for speeding solution 
by successive approximations may be used to reach finally a set of values 
which satisfies all the equations. 

Without doing this calculation, we may predict that it is unwise to use 
values of /,, and /,, which differ too much, since the heat-transfer rates from 
the two portions must be equal. 

Here the calculations will be done only for the case when both portions 
of the tube have Grashof numbers above the critical, so that N, oc G}, 
N, x (l/a)Gt. Then, by (93), 





Gac _ (Sefer ere sore |t (94) 
: = ae ; ( 
Gan Saka Yn Kala 
and it follows by (91) that 
q 7 > 2 78 y.. yee 
Gan = Tu —T, nae 14 Gas YH <u!) (95) 
’ ‘ 7 , ae 5 fi 
Gon Tn 0c Lo So ké lé:/ve Ke 


(As suggested in section 1, this is } when all constants are uniform.) 


Hence the Nusselt number based on the total temperature difference, 


Y Y 
NC Q ae Gan lt 9.) 1G Gan 
Vin anak, 7 TT.) YH Ge 7aH GC 
27k y (Lo — Loc) aH a 7aH 





0-1] ln (Go) ad (ett Vy <u! 3 (96) 
e*”* ; 


23 73 |, 
ro fc kb ld/ve ke 


If all constants are given, including /,,+-/, but not the ratio /,,/1, then (96) is 
a maximum if 





= (222-2 C (97) 
lo Me Sa ki /Yn Kn 


ln ( roto ke: tet) 
which will not be very far from 1 in practice. 

There will also be a condition on /;,//, in order that G,,, and G,, shall be 
supercritical with as much as possible to spare. This means that the smaller 
of G47 (lq7/a)-™ and G,,¢(/,/a)-™ shall be as large as possible, where M is the 
slope of log(G@, .»it) against log(//a), for values of l//a near half the total 
length-radius ratio; thus MW is three less than the slope of the curve of 
Fig. 10 at this point. Now it is fairly clear that the two quantities alter in 
opposite directions as /,,//, is altered, so that the maximum minimum is 
achieved when they are equal, which gives 


bi 3(4M +38) 
ln = (‘at key Up <u) , 


lo ao fol 





(98) 
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Since M lies between | and 8, (98) also means a value of near 1, especially 
for large total lengths. Thus although (97) and (98) would normally be 
incompatible, one could satisfy the more important one, namely (98), with- 
out departing very far from the maximum heat transfer condition given by 
the other. 
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NOTE ADDED IN PROOF 


Si>ce the above was written Mr. B. W. Martin, Lecturer in Mechanical Engineering 
at King’s College, Newcastle, has done some experiments on the system studied in 
sections 2-12, with glycerine as the working fluid and l/a 47-5. The modified 
Grashof number G,, based on tube radius and on the values of the fluid properties 
at the temperature of the wall, ranged from about 600 to 3.10° Turbulence is 


believed to have begun to appear only at the highest values of G, that were attained, 
and to have been absent at G 1-5.10°. This seems to indicate that transition 


values of G, are a lot higher than has been assumed in the paper, although it should 
be noticed that these values of G, are based on the viscosity of the very hot glycerine 
at the wall, and that the cold glycerine at the centre of the tube had a much higher 
viscosity, which will have tended to damp turbulence starting up there. 

The lower curve in Fig. 7 is confirmed in shape for 1-0 logigt,; < 3-4, but the 
experimental points lie above the curve by about 0-25. The upper curve is well 
confirmed for 4-3 logit; 1-7. Between these ranges the flow exhibited regular 
periodic oscillations which are believed to be associated with the existence of the 
two possible steady flow régimes (according as the boundary layer does or does nct 
fill the tube) of sections 4 and 6; these oscillations sometimes occurred for logy, ¢, 
as low as 3. The amplitude was greatest around logy tf, 3-5, beyond which it 
decreased as ¢, increased, until at logy,¢, 4-3 (that is, at G, 10°) the oscillations 
had disappeared, although it is believed that turbulence had not yet set in. 

Thus the presence of two curves in Fig. 7 appears to represent not a mathematical 
error but a genuine physical uncertainty, which is resolved in practice by the 
appearance of an oscillatory flow. However, it remains to be seen whether the above 
interpretation of Mr. Martin’s experiments is confirmed by further study. His 
investigations will be published when they have reached a more advanced and 
conclusive stage; the author is grateful for being permitted to include the above 


brief reference to them. 





THE PROPAGATION OF SHOCK WAVES IN A 
CHANNEL OF NON-UNIFORM WIDTH 
By W. CHESTER (Dept. of Mathematics, The University, Bristol) 
[Received 17 July 1952] 


SUMMARY 

The disturbance produced behind a plane shock wave of arbitrary strength, 
travelling down a channel of non-uniform width, is investigated. The problem is 
linearized on the basis of small variations in the width of the channel. 

It is found that, if the initial and final widths of the channel are 26, 2c, respectively, 
the ultimate effect in the neighbourhood of the incident shock is simply a change in 
the shock strength with a consequent change in pressure behind the shock of arnount 

K(p,—po)(1—e/b), 


where (p; — Pp) is the initial pressure discontinuity across the shock. The parameter K 


decreases monotonically with the shock strength and 0-5 > K > 0-394. 
When the flow behind the shock is subsonic an acoustic pulse also travels back 
along the channel with sonic velocity relative to the fluid. If 2g(a) (= 0, if 2 < 0) is 


the variation in the width of the channel, the pressure in this pulse is ultimately 
Nb-(p,— po)g(78), 

where s is the distance from the front of the pulse. The parameter 7 decreases mono- 

tonically with the shock strength and 2 > 7 > 0-608; N takes the value 0-5 for weak 

shocks and is singular when the flow behind the shock is sonic. 

When the flow behind the shock is supersonic a similar pulse is produced (N becomes 
negative) now travelling downstream, and this is superimposed on the solution for 
steady flow in the channel. 

The singularity in N for sonic flow behind the incident shock would suggest the 
possibility of strong secondary shocks occurring. 


1. Introduction 

THE object of this paper is to investigate the flow which results when a 
shock wave travels along a two-dimensional channel of non-uniform width. 
The variations in width are assumed to take place within a finite length of 
the channel, and this region of transition separates two portions of infinite 
length and uniform, but not necessarily equal, widths. The problem is 
linearized on the basis of small variations in the channel width, but there 
is no restriction on the shock strength. The pressure field behind the shock 
front is built up from the known solution of the diffraction of a shock wave 
travelling along a wall when the direction of the latter changes discon- 
tinuously, so presenting a corner to the flow (see Fig. 1 and (1)). 
From this is deduced a basic solution for the flow in a channel formed by 
two such walls, that is, a channel of uniform width connected to a channel 
of uniformly increasing width (see Fig. 2). Because of the linear character 


of the solution, the flow in the case when the cross-section varies in an 


{Quart. Journ. Mech. and Applied Math., Vol. VI, Pt. 4 (1953)] 
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THE PROPAGATION OF SHOCK WAVES 441 
arbitrary manner can be obtained by superposition. This solution is then 
used to investigate the ultimate nature of the flow at large distances from 


the transition section. 


2. Diffraction by a corner 

We consider first the two-dimensional problem envisaged in Fig. 1. A 
plane shock moves with uniform velocity U into still air and is bounded bya 
plane wall perpendicular to the shock front. A second plane wall meets the 


incident 
ly ~ shock 








/ 














Fic. 1 (6) 


first along a line parallel to the shock front, the resulting corner presenting 
an angle (7+-8) to the flow. When the shock front has progressed beyond 
the corner it will be diffracted by the inclined wall, and an expanding 
region of non-uniform flow will be formed behind the shock front. 

Let the fluid velocity, pressure, density, and sonic velocity in the region 
of uniform flow behind the shock front be denoted by q,, , p,, @, respec- 
tively, and let the corresponding quantities ahead of the shock be 0, po, po, do. 
Choose coordinate axes O(x, y) with origin at the corner and the z-axis along 
the original wall produced, and let the shock reach the corner at time t = 0. 
If gq, < a,, then at any later time the perturbed flow is bounded by the 
incident shock and a circle of radius a, ¢ with its centre at the point (q,¢, 0). 


Ifthe uniform flow behind the shock is supersonic, q, > a@,, then the corner 
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lies outside the sonic circle and the boundary must be augmented by the 
tangent to the circle from the corner. In either case, since there is no 
fundamental length in the data defining the problem, the flow variables are 
expressible as functions of 2/t, y/t only. 

If the angle 6 is small, the pressure in the region of perturbed flow will 
differ from p, by a term proportional to 5 on a linear theory; this excess 
pressure is conveniently denoted by 


(p,- pao" an! de 5 = (p,— Pp) p(X, Y) 5 (say). (1) 


at a,t 


The linearized problem has been solved by Lighthill (1), and the following 
is a summary of his results. 

Let M (= U/a, > 1) be the Mach number of the shock, and let ¥, 
(= q,/a,) denote the Mach number of the uniform flow behind the shock. 
Then 

M, = -= Fe) one ] (2) 
| (7M?—1)(M?-+5)]! 

(L. (3)t; the adiabatic index is given the value 1-4 appropriate to a 
diatomic gas). 

In the (X, Y)-plane the circular part of the boundary has unit radius and 
centre at the origin. The corner is at the point (—,, 0), and the straight 
part of the incident shock is the line 


U—q _ ( M?-+ 5 \3 


A = = * = COS (3 
ay \7M2—1) x ) 
(L. (9)). 
The excess pressure coefficient p satisfies the relation 
Cp, op _ C{D(z,—1 a tS (4) 
CY, . CX, (zi 1)*(z,—1 T Pa 1) [B—i(z, 1)] 


(L. (45)), where 


l/ (pexpi(O+x)}—1\*_, [ pexpt#—expiy |? . 
24 ety, = — » (5) 
2|| pexpi#—expiy | | \pexp{i(@+ x)}—1) | 
ae. Se (6) 
1+ p? 


+ This denotes equation (3) of Lighthill’s paper (1). To facilitate comparison, Lighthill’s 
notation has been used in the main, but the following changes should be noted: 
(x, y), (X, Y), k and p now become (X, Y), (x, y), cos x and 8p cos x, respectively. 
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(L. p. 459 et seq. and (31), (35)), 


B= 2M, atP=2!M?cosx,  y* = 2(M,+cos x)*(Af, cos x+ 1), 
(7) 
0 2M,(.M/,+-cos x)(B ty My T x)sin x (8) 
a(M, cos x+-1)? : 
5(.M, cos x+- 1)? cos x D(y+a)(a4 B) _ x tB+y (9) 


6.MF(M, COS x) xB xBy 
(L. (3), (9), (15), (50), (51).) 


3. Application to flow in a channel 
We now consider the propagation of a shock between two walls of infinite 
extent forming a channel of uniform width 2b (say) connected to a channel 


of uniformly increasing width (Fig. 2). In the initial stages the solution of 











Fic. 2. 


section 2 is still valid, the diffraction pattern produced by a single corner 
being repeated at the upper and lower walls of the channel. This will be 
the complete solution until the circular waves have expanded sufficiently 
to encounter the opposite walls, when further modification is required to 
satisfy the boundary condition there. If, for the moment, we ignore the 
divergence of the channel walls (in so far as the reflecting property of the 
latter is concerned) and imagine them to be straight and parallel along 
the whole length of the channel, then the modification would be a simple 
reflection at each wall simulated by appropriate images. This is also 
correct to the first order if the divergence of the walls is small, for the flow 
quantities in the perturbation waves are already small quantities of the 
first order. 

The image waves will themselves be reflected eventually; thus the com- 
plete linear solution is obtained from an infinite image series expanding in 
time with sonic velocity, their centres moving downstream with the velocity 
of the fluid behind the original shock. 

Let the w-axis now be taken along the centre-line of the channel, and let 
the y-axis separate the two sections of uniform and non-uniform width. 
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In terms of the basic solution of section 2, the excess pressure inside the 
channel is then given by 
" = = (x—qyt (2n+1)b—y) | 
8(p,—p,) P(x, y,t) = 8 > lyin ec As 
(Pi— Po) F(x, y, t) = 0(P1— Po) P| a,t a,t | 


n 
seth! o+y)\ (10) 
| a,t a,t | 
where the summation is taken over all the images which have expanded 
sufficiently to influence the pressure at the point considered. 

To generalize this for an arbitrary variation in channel! width, we first 
consider the effect of moving the origin of coordinates upstream a distance €, 
so that the channel has constant width for —o < x < &. If the origin of 
time is suitably altered so that the shock reaches x = 0 at t = 0, the excess 
pressure inside the channel will be 


5(p,;— po) P(a—€, y, t—E/U). (11) 
Now let the upper wall of the channel be defined by 
y = b+ g(x) (12) 


(the channel itself being symmetrical about the centre line), and choose the 
origin so that all variations in width occur for x > 0, so that g(a) = 0 for 
<= U. 


Equation (12) is equivalent to 


dy/dx - ( H(x—€&) dg'(é) = g'(x), (13) 
0 
where 
H(x l (x>0 
(x) (a ) (14) 


and is the Heaviside unit function. 
Since the expression (11) is the excess pressure appropriate to a channel 


for which dy|dx = 8H(x—8), (15) 


we identify dg’(é) of (13) with 6 of (15), and sum for all the elementary 
variations in slope to get the following expression for the excess pressure, 
(Pi—Ppo) | P(x- E,y,t g U) dg'(&). (16) 
0 

4. Asymptotic behaviour near the shock front 
The main purpose of this paper is a study of the asymptotic behaviour 
of (16) at large distances from the origin when the variations in the width 
of the channel occur within a transition section of finite length joining two 
sections of uniform (though not necessarily equal) width. 
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Equation (10) is first transformed by means of Dirichlet’s summation 
formula, namely 
ae, 
>” f(r) lim > | f (z)e?™*" dz, (17) 
Vi n Ve Anu n —Ay 
where =” denotes that the first term is }f(v,) if v, is an integer, and the last 
term 4}/f(v,) if v, is an integer (3). 
We then have 


P Lp wv nt 6b i) 11 (= nt 6 ry\ 
1,¢ a, 


- a,t ( t a,t 
1 v—qit)/ayt . 
a,t ~ . 1) a Vlex {mare ‘Y 4 \ay- 
tp J Oe ee 
? (b—y)/ayt 
l r—qit)/ayl 
. fjr—q,t ,) (nae : ) oy 
1)” Vie tY — 1Y}. (18 
), Oat ee 
(b+y)/ayt 


We consider the behaviour of (18) at a fixed distance behind the main 
shock, so that 


x Ut—constant. (19) 
Then, when ¢ is large. 
Sin xX ° r 
Pr = fi (—1)"pjcos x, Yiexpo™ cu cos ""4 d¥+-A4 o(;) 
b Land B b ) t 
giants (20) 


where A is some constant independent of t. 

Comparison with the result for diffraction by a single corner (section 2) 
shows that the integral occurring in (20) is taken along the shock from the 
wall to the sonic circle. (This follows from (3) and the fact that in the (X, Y)- 
plane the sonic circle has unit radius.) At the upper limit p = 0, the pressure 
being continuous at the sonic circle to the present approximation. Equation 
(20) thus becomes, after an integration by parts, 


n x sin x 


/ t ‘ 4 n oa fe : 1 
P “I Y dp—2 a ae cos 74 | sin” ™* dp+A-+O(t-). 
b Lan NT b b 
) n l 0 
(21) 
Now on the shock Yy 0. 
x 1\3 . 
} (= }* sin x (22) 
%,+1, 


(equations (5) and (6) or L. (47)), thus x, runs from 1 to oo. Further, from 
(4), ‘f a 
Op C| D(x,—1+y?) 1 |(a +-B)sec y 


y 
x, (a, +1)(a,—14-y?)(a,—1+02)(x, — 1+?) 


(23) 
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Thus, since the second term in (21) is O(t-"), it follows that 
: — st ; 
sin x(2 I) cP dx,+A+O 
X,+1) oz, t 


" , D(a? +-y?)—1 |x? dx 
— 20a, Na+ Bytan x | a5 He ye + aS ty?) 
é 


+ A+ o(;) 
t 


—KUt/b+A+ O(t-*) (say), (24) 





- —a,tb-} 


« 





where 

K — (Cr(a+B)tan x}{(Dy?—1)(2!+a+B+-y)+ D(2toB+ hay +2'By + aBy)} 
(M, + cos x)(2!+-a)(2!+8)(2!+y)(a+B)(a+y)(B+y) 

(25) 

Finally, when the variation in cross-section is arbitrary, the pressure 

change behind the shock wave is obtained by substituting (24) in (16). The 

range of integration is effectively finite in (16) because of the assumption 





that g’(2) is identically zero outside a finite interval in the neighbourhood 
of the origin. The pressure change is thus 





—(r—pa) { [A t-g/0)-A+ 00 | dg'(é), (26) 
) 
and, since g’(€) is zero at the end points, this is 
: —_ Po) rg) | Olt 1), (27 
y) 


where |g] denotes the net change in g. If the final width of the channel is 2c, 
then (18) may be written as 
K(p,—p»)(1—eb-!)4+- O(t-?). (28) 

Thus, ultimately, the effect of the change in cross-section amounts to a 
change in the shock strength. The flow behind the shock again becomes 
uniform, the pressure (in excess of p,) being given by (28). 

The parameter AK is a monotonic decreasing function of the shock 
strength. The limiting values for a shock of zero strength and a shock of 
infinite strength are 0-5 and 0-394... respectively. The detailed variation 
of A is shown in Fig. 3, where it is graphed against p,/pp. 





A simple verification of the result for the limiting case of an acoustic 
pulse (K = $) is possible (2). For the energy in a progressive wave, per 
unit cross-section of the channel, is proportional to the square of the 
amplitude. Since the amplitude of the reflected pulse will be O(1—eb") | 
it follows that the energy in the transmitted pulse will be unchanged to 
the first order in (1—cb-"), and this will be the case if the amplitude varies | 
inversely as the square root of the width of the channel. The excess pressure 
in the pulse is a convenient measure of the amplitude, and initially this is 
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(py—Po). The ultimate value is thus 


(Pp, P_)(b/c)? ();—Po)-+t 4(p, Po)(1—e b)+O(1—ce/b)?. 
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Fic. 3. 


Asymptotic behaviour of the returning pulse 

When the flow behind the incident shock is subsonic, a disturbance will 
also be propagated into the upstream section of the channel. At a large 
distance from the origin the wave-fronts of the original diffracted waves 
and their subsequent reflections will have an envelope consisting of a plane 
wavefront travelling with sonic velocity relative to the fluid. At any 
given time the position of the envelope will be given by 

r (a,—q,)t, (29) 

and behind it the flow will tend to become uniform across the width of the 
channel. 

We find the pressure variation in this returning wave by considering the 
behaviour of (18) at a fixed distance (d, say) downstream of the envelope 
which constitutes the front of the wave, so that 


v (a, q,)tA d. (30) 
Then 
P sp d b ) , 7 1 d b4 ) | 
a,t a,f a,t at, 
- {1—(1—djayt ; : 
yt wo : = ( \ (nm i ay 
* 2. | | { x0 9} P|- - at” Y exp} , (a, t) ") 1) 
1—(1—d/ayt)*} ,; 
( nr - 
1)"p{ 14", Y tY—y)ldy|. (31 
| (—1) pI at Jexp| , Wie (31) 
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We therefore require the behaviour of p(X, Y) when 

X =m —t4-OH-, 

Y = O(t-'). (32) 
It may be deduced from (5) and (6) that, when X and Y satisfy these 
conditions, 


a = —14+4{X+1—Y741+4 O(t-!)}tan?y/2 | e 
y, = 4¥{(2X+2—Y°H14 O(t-)}tan2y/2 J (8) 
Further, in the neighbourhood of z, = —1, 
p — —OPA=*)+ +a tte Pty P+ Of + ay) 9) 
cos y(2—+y?)(a-+2*)(B+ 2!) 
(34) 
This follows at once from (4) and the fact that p = 0 at z, —1 since 


this point is on the boundary of the diffracted pulse. A combination of 
(33) and (34) gives 
p = —L[1+4+ O(t-!)|[2(X¥+1)—Y?} (35) 
_ 20[D(2—y?)+ 1]tan x/2 
cos x(2—y*)(a+ 2!)(B+2#)’ 
when X and ¥ satisfy (32). 
Equation (35) for p is now substituted in (31) to give 


La,t 
? 4+ O(t-3)]s. 
I P [1+ O(t-*)| 


with L (36) 


(2d/a,t)* 
- 


m / 9 1 , 4 = 
S (—i7 ~~ Y2)"exp = cos wu) ay + O(t-*) 
} a,t b b 

0 


Lrd/2b+-O(t-*) 
— Ln{x+(a,—q,)t}(2b)-1+ O(t-4) (37) 
in virtue of equation (30). 

Finally we substitute (37) in (16) to find the pressure variation in the 
returning pulse when there is an arbitrary variation in the cross-section 
of the channel. Since P is clearly zero ahead of the returning pulse 
({v+-(a,—q,)t} < 0) it follows that (16) becomes 


oh (P1—Ppo) | fy—£+(a,—q,)(t—E/U)} dg’ (€)+- O(t-') 


0 


U{ay+ (ay—q)t}/(U +a,—q1) 


(Nb-!)(p,—p,)gl rit (a,—9q,)}|4+O0(t), (38) 
U _ M,+cos x 


where T= — _= — 
a,—qt+l 1+-cos y 


(39) 
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ry 
Tr ( 





D(2--y*)+1]tan x — (40) 


( N Li /2r a 
and (2—y?)(a+2#)(B+ 2#)(.M, + cos x) 


Thus a graphical picture of the ultimate nature of the pressure variations 
in this rearward pulse is obtained by imagining the upper wall of the channel 
stretched in the x-direction by a factor 7, magnified in the y-direction by 
a factor —N(p,—p )/b and moving upstream with a velocity (a,;—q,). 















































4. |} —___—_- 


The values of 7 and N for various shock strengths are shown in Figs. 4 
and 5. The parameter 7 is monotonic decreasing and takes the values 2, 1, 
and 0-6076... respectively when the incident shock has zero strength, when 
the flow velocity behind the shock is sonic, and when the shock strength 
becomes infinite. 

The parameter NV is not so well behaved. For subsonic flow behind the 
incident shock it increases rapidly from a value of 0-5 for a shock of zero 
strength, and becomes singular when the flow behind the shock is sonic 
Near M, | the behaviour of .V is, in fact, described by the asymptotic 
relation | 

N ~ . (41) 
(1—.W})cos x 

The solution, of course, breaks down near M, 1. However, experience 
of other linearized solutions suggests that large non-linear disturbances 
are often concealed by such singularities. If this is a correct inference in 
the present case, then a strong secondary shock wave ought to be formed 
if there is a contraction in the channel and the flow behind the incident 
shock is just subsonic. Thus we obtain a mathematical expression for the 
process of ‘choking’. 

24 


5092 
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For supersonic flow behind the shock N becomes negative (and singular, 
according to (41), near MV, 1). As the shock strength tends to infinity 
the value of N tends to 1-550.... 















































In the supersonic case it must be remembered that although (38) can 
still be interpreted as an acoustic pulse (now travelling downstream with 
velocity (q¢,—4a,)), it no longer represents the whole of the disturbance. 
In the case of a uniformly expanding channel there is the additional effect 
of the linear approximation to the two Prandtl—Meyer expansions at the 
corners (see Fig. 1 (b)), and these are themselves continually reflected by 
the channel walls. In the general case when the variation in the channel 
width is arbitrary, the additional disturbance can be thought of as the 
integrated effect of these Prantdl|—Meyer expansions (or weak compression 
waves) and their reflections, which originate at each elementary incre- 
ment of slope in the walls of the channel. Evidently the mathematical 
solution representing the combination of all these must be simply that of 
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steady flow along the channel, and if the velocity vector is defined in terms 


of a velocity potential ¢ as q,{1 ns Po A then ¢ satisfies the equation 
Ga’ ey) 
24 r n2 
P — (M?—1)< ¢ (42) 
oy” ox" 
a i 
and a boundary condition m2 “4 (43) 


cy dx 
on b 
In the Heaviside operational calculus, (42) and (43) become (with s as 


the operatol! ) 


C-@O 


2 — (M2—1)s’¢6 = Bs*d (44) 
Ccy- 
ad F de 
say), and s | e-F dx » (say) on y = b. (45) 
cy - dx 
The solution, symmetrical about y = 0, is 
d Bs)cosh Bsy sech Bsb. (46) 
\lso, by Bernoulli’s theorem, the excess pressure is, to the first order 
4} 
A = —p, gi = —pigis¢, (47) 
OX 


and so, by interpreting (47), we have 


A Pit : nan _ f sx ds, (48) 
27i ~=})~=6Bs sinh Bsb 


Now dg/dx is zero except in the neighbourhood of the origin. Thus the upper 


limit in the integral occurring in (45) is effectively finite and, for small s, 


|g |s+-O(s?). (49) 
It follows that, if » is not singular to the left of the contour of integration 
in (48). : 
dz “ nat NrX 
A Pitiry, S A, cos| “Neos (50) 
Bb ya b | B 
7 
for some constants A,. The series in (50) represents a fluctuating con- 


tribution to the pressure, and this fluctuation does not die out even at large 
distances from the origin of the disturbance. 
It can easily be shown (with the help of L. (1)) that, in terms of the 


pressure discontinuity across the incident shock, the first term in (50) can 
be written 
M,(~;—Po) (g] 
(M? 


1 


(51) 


1)b cos x~ 
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and this, in combination with (38), represents the ultimate average pressure 
increment in the pulse. This term is also singular at MV, 1, and com- 
parison with (38) and (41) shows that it is asymptotically equal and opposite 
to the value taken by (38) in the heart of the rearward pulse (when 
g|r{x-+-(a,—q,)t}] takes its final and constant value, here denoted by [g}). 
Thus, in combination, the singular terms tend to cancel out in the heart of 
the wave, but the linear theory still indicates rapid increase of pressure 
(if [g] > 0) as the pulse is penetrated from its front. Hence, for reasons 
similar to those invoked in the subsonic case, there should be a secondary 
shock produced when there is an expansion in the channel and the flow is 
just supersonic behind the incident shock (secondary shocks in expanding 
channels have, in fact, been observed—see (4)). 

It should be mentioned that the above reasoning does not apply to the 
flow immediately behind the incident shock, where the steady flow field 
has not penetrated. Equation (28) is carried over to the supersonic case 
without qualification. 
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ON SOME EIGENVALUE PROBLEMS OF EXCEPTIONAL 
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SUMMARY 


Relaxation methods were employed in 1941 (2) to estimate the critical load 
ind the mode of distortion (‘waving’) for a flat plate representative of a cantilever 
I-beam. The ilue thus obtained for the critical load was believed to be correct 
within 1 per cent.; but the mode was less accurately determined, and attempts to 
mprove it have resulted in the present paper. (Though less important, practically, 
than the critical load, it has interest as determining the arrangement of stiffeners.) 

Normally g. for struts or ‘whirling’ shafts) the eigenvalues are widely separated 
und all positin and the gravest mode can be guessed with fairly close accuracy 
since it is devoid of nodal lines. But in general these features cannot be postulated, 
ind then approximate treatment is much harder. All three appear to be lacking in 
he cantilever problem, which thus exemplifies problems of a wider class. The conse- 
juent difficulties, being mainly computational, are not met by discussions which 
presume that equations can be solved exactly. 

Introduction 


1. THE complete solution of an eigenvalue problem (giving every eigenvalue 
and its associated mode) is seldom obtainable and from a practical stand- 
point is seldom needed. In a case of elastic instability, for example, what 
matter are the load which first induces instability and the mode of distortion 

waving’) which results. The first is usually estimated on the basis of 
‘Rayleigh’s principle’, the second by a process (here termed ‘intensification’ ) 
which is commonly attributed to Vianello and/or Stodola.+ 

Applied to systems like struts or ‘whirling’ shafts, these methods quickly 
yield the wanted information: Rayleigh’s, because the deflexion in the 
gravest mode is known to be one-signed and as such can be guessed closely; 
‘Intensification’, because the eigenvalues (critical loads or speeds) all have 
the same sign and are, moreover, widely separated. But in general these 


conditions are not realized, and on that account both methods (as normally 


ne questions of attribution are discussed in section 34. 


[Quart. Journ. Mech. and Applied Math., Vol. VI, Pt. 4 (1953)] 
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applied) may fail. Thus the few solved cases of elastic instability in plates 
suffice to show that in such problemst 

(A) the gravest mode may be characterized by nodal lines, and as 

such be impossible to guess with any accuracy; 

(B) the eigenvalues may have both signs; 

(C) two or more of the gravest eigenvalues may be nearly or exactly 

equal. 
Any or all of (A)-(C) may be confronted in cases which, on account of their 
complexity, are tractable only by numerical methods; and must be sur- 
mounted before a relaxational technique can be applied in the manner (e.g.) 
of Pellew and Southwell (1). This is the object of the present paper. 

2. We take as a test example the system shown in Fig. 1—a problem 
attacked by L. Fox, J. R. Green, and F. 8. Shaw in 1941 (2). A deep 
rectangular cantilever, with edges ‘clamped’, is stressed by uniform thrust 
applied to one of its longer sides. Contours in Fig. 1 exhibit the compressive 
principal stress which this load induces when acting downwards: then, the 
gravest mode is characterized (Fig. 2) by waving localized in the bottom 
right-hand corner; but clearly, if the loading were reversed, the principal 
stress would be compressive, and the waving localized, in the upper right- 
hand corner. The system thus exemplifies (B), section 1. 

The equation defining the characteristic modes and eigenvalues can be 
deduced in the manner of (3) sections 494-7. It has the form 

ow AD’ w, (2)* (1) 
where 


ow DV4w (v: te ! = | 


: C ow (OW c Cow OW 
Ad w : (P, — ; = (P, 7 a > 
OX \ ~ ox cy cy\ ~ cy Cx 


\ 


(2) 


(D denoting the flexural rigidity). For the loading (w) shown in Fig. 1 
(cf. (3), section 416) 


Pp: al 32°(2y d)—4y3 + 6dy? —12d?y-+-1d3 | 

Fy 
) W (998 2da2 9 
P, a (2Y" 3dy*) big (3) 
7 

6 ‘ 
S = 5 xy(d—y) 


+ A and C are exemplified by a rectangular plate compressed in one direction; B by a 
square or rectangular plate which ‘buckles’ under shear. Perhaps because only simple cases 
have been solved, it seems not to have been noticed that in general a plate can buckle 


under the load-system either as given or reversed. The theoretical significance of this fact 
was first remarked by D. N. de G. Allen. 
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(Contours show the (compressive) principal stress as a multiple of w/2h.) 
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; wd? 
so the eigenvalue Ax Ee (4) 
Dh 


when (1) is made ‘non-dimensional’ by the substitutions 


x =: dz’, y = dy’. 


3. Weshallsuppose (1) to have been thus rendered ‘non-dimensional’ ,and 
in treating the test example we shall use 3, 3’, indiscriminately, both for the 
operators specified in (2) and for their approximations in finite differences 
((2), sections 24-25). But the symbols admit of a wider interpretation, and 
thereby the range of this paper is extended:+ they may stand for any two 
self-adjoint operators, provided that in the ‘energy equation’ derived 
from (1)—viz. 

I (say) | wdwdady = r || wi'wdady = Al'/A (say), (5) 

the integral (I) on the left is necessarily positive. This can be shown to be 
the fact when (as here) # = V4 and the edge-conditions are 

ow : 

w 0. (6) 


Cv 


No restriction is imposed on the right-hand integral (/'/A) in (5), so 
the ergenvalue X may have either sign. 


Assumptions 
4. Asalways (when itis linear) the governing equation does not define the 
amplitudes of the characteristic modes. Typifying them by w,,, a function 
defined by , ‘ 
: Pw, Xi. o Ww}. (7) 
and by the conditions (6), we ‘normalize’ them by imposing the relation 
|| w,, Pw, dady Ay | | wy, Fwy, dady=1, (3)* (8) 

and we presume them to satisfy the ‘conjugate relations’ 

| | w, dw, dady | | w, dw, dady = 0] 


ae - . (4)* (9) 
| | wo’ w, dady | | w, Ww, dady of 


Finally, we assume the validity of expansion in the characteristic modes: i.e. 


that, when the A’s are given appropriate values, the series expression 


w > [ A, wy; | (5)* (10) 
r 


Its application depends, of course, on whether (1) is tractable relaxationally. 
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holds in respect of any mode which does not violate the edge-conditions, 











4) and admits of operation either by & or by #. 

(Only this last assumption is open to question: strictly, it calls for 
examination in respect of each particular equation. But it is justified in 
relation to systems of finite freedom, and in effect we change the system 
to one of this kind when we substitute finite differences for derivatives in 

id the governing equation, thereby reducing the number of the points that 
he are ‘relaxable’. This standpoint was adopted in a recent paper (4) which 
os gives proofs of some theorems used here (e.g. that the assumption J << 0 
id requires all the eigenvalues to be real). Those theorems will be cited without 
si proof but (to facilitate reference) with their equation and section numbers 
d in (4) distinguished by asterisks.) 
5) I. REVIEW OF APPROXIMATE METHODS ALREADY CURRENT 
™ 5. Taking the test example as representative of problems which may have 
| any or all of the features ‘A’, ‘B’, and ‘C’, section 1, we now examine the 
applicability of some current methods. No hope of exact solution can be 
6) entertained, so we restrict the review to methods which in effect replace a 
continuous system by one having finite freedom. Equally it is evident that 
SO a mode like Fig. 2 (Fox and Green’s solution (2) of the test example) will 
not be adequately represented—as in methods like Galerkin’s—by a 
synthesis of a few selected functions chosen without prior knowledge of its 
| jorm. In Fig. 2 the number of nodes that are ‘relaxable’, and hence the 
he | number of the characteristic modes and eigenvalues, is 273. It appears 
- that little has been written with systems of this complexity in mind; and 
perhaps for that reason little attention has been given to the consequences 
) of inescapable computational errors (cf. (4), section 7*). 
. Rayleigh’s method 
(8) 6. On the assumptions made in section 4 
[ || witw dady > (A?) | 
hi (6)* (11) 
I’ /A ( [ wu Ww dady > (Az A;.) 
) Js 
| consequently (section 3*) Ap as defined by 
e. | Ap = All|! ¥ (A2)/¥ (. Z/A,.) (7)* (12) 
has a stationary value when the mode w is ‘normal’ so that all of the A’s 
0) | are zero except one. This is ‘Rayleigh’s principle’—propounded by him as 


holding when the eigenvalues are all positive. It holds also when (as here) 


they have both signs but (section 2*) are all real. 


R. 
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(Absolute magnitude of distortion is indeterminate.) 
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Lord Rayleigh (5) argued, further, that Ap as defined in (12) always has 
a value intermediate between the least and the greatest of them (or above 
the least, if the system has infinite freedom), consequently Ap as computed 
from a guessed form w (here termed a ‘Rayleigh estimate’) may be regarded 
as an upper limit to A,. This theorem fails when /’/A can take either sign 
or vanish, for in the last event A, is infinite. But a theorem due to G. V. 
(section 4*) asserts that 

Ap, cannot lie between the least positive and the least 
negative eigenvalue (A, and A_,, say).t (13) 
‘Intensification’ 

7. The practical importance of Rayleigh’s principle lies in his use of it 
section 6) to obtain an upper limit to A, by simple computations based on 
a guess regarding the form of w,. In most problems of the sort which he 
considered the form can be guessed with considerable accuracy, and then 
in virtue of the stationary property the resulting estimate of A, is extremely 
close. But in the wider class considered here little confidence can be put 
in a guess regarding w,: then it may be improved by a process, here termed 
intensification’, which is variously ascribed to Schwarz (6), Vianello (7), 
ind Stodola (8). This derives from w, successively, other forms wy, wyy,.... 
ete., satisfying the edge-conditions together with 

dw, = Ad'w | 


VW yy A ow 


in which A is any convenient multiplier. 


(9)* (14) 


Provided that w, is a component of the guessed form w, when every eigen- 
value > 0 then wy, wyy,..., etc., will approximate more and more closely 
to w,; when the eigenvalues have both signs they will approximate to w, 


or w_, according as A, or A_, is the smaller. For if w is expressible in the 


series as s (A, Wy), (10)bis 
then Wy VS (Aj, w,/Ax) | 
Wry \2 S (A; w,/AZ) (10)* (15) 


, etc. 


so the gravest mode becomes successively more predominant in wy, Wyy,..., 


ete. 

\ normal mode will be termed positive or negative in conformity with the sign of its 
associated \. © will denote a summation extending only to positive modes, typified by w, ; 
x a summation extending only to negative modes, typified by w_;; X a summation 


ending to every normal mode 
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8. Another interpretation is admissible. The response (w) to a transverse 
loading Z, is governed by 
[I—A # |w = Z. = Huy (say), (16) 


when the edge-loading has an intensity defined by A; and a formal solution 


of (16) is w > | A, w,./(1 A A,)}) 


ra ’ (17) 
when Wo 2 | A; W,.| 
Alternatively, assuming w to be representable by the series 
w= Pot+PhAt+...+P,A"+..., ete., (18) 
we obtain > [A(aP, —# P,,_,)| = dw, 
on substitution in (16); and hence 
oP, oP. bs OP, Fwy. (19) 


From these last relations and from the edge-conditions it follows that 
P= w, = > [A,w,}, P > [ A, w;,/Az], 
so that when » is large (and if every A, \1)s 
P, > A,w,/A} (therefore A, P, > P,_,) as A> A,. 

Then the series (18) > A, w,/(1—A/A,)—in conformity with (17). 

Thus successive intensifications of w, yield cofactors of A, A?,..., etc., in 
a power-series (18) valid when A < A,, and as such they have the nature 
of successive corrections. (In section 7 they appeared as successive approxi- 


mations.) 


‘Partial transfer’ 
9. That A, < |A_,| in the test example can be shown by resolving the 
load-system of Fig. 1 into two components, (a) and (b), Fig. 3. If (a) were 


tw sw 
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to act alone, the eigenvalues would form equal and opposite pairs; but 
superposition of (b), a destabilizing system, will alter each of them by a 


positive quantity, thus increasing A,| and decreasing |A_,). 
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erse Here, then, ‘intensification’ may be expected to furnish an estimate of 
the wanted mode w,. But ‘B’, section 1, will frustrate that purpose if 


\ 


16 r_, A,, so may be a cause of difficulty. In that event a device which 
\ ) 


we here term ‘partial transfer’ may have value, as giving a bias towards ; 
“ion components of the wanted sign. (An equivalent device was proposed by 
L. F. Richardson (9) in 1910.) Its success requires the proviso of section 7. 
(17 A, denoting some chosen number, and A’ standing for (A—A,), then (1) 
may be written in the equivalent form 
19 Xo P |w rN Pw. (1)A 
8 
Let the operator on the left of (1)A be substituted for # in (14), so that ‘ 
Wy}, Wyz,---, ete., are now related with w by 
I—Ay H |wy A Fw) 
: (20) 
(19 vey eC. 
Then an argument on the lines of section 7 shows that 
wy \ > | A. w,,/ (Aj, Xo) | 
Wry ass | A, Wy (A;. a (21) 
» etc. 
when w has the series expression > | A;,w,|: ie. A, in (15) is now replaced 
_ aa . . ° ° 
by (A;.—Ap), So now the measure of the relative contributions of w, and w, 
fure th ‘ ; = . . . . 
is reduced in each intensification by a factor (A,;—Ag)/(A;,—Ap,) instead of 
Or) ? on ° 
A,/A;, and convergence will be much accelerated if A, = A,. In particular 
negative components will be eliminated, because |(A,—A,)| > |A;,| when 
1. < 0. 
the iia ‘ ; : ‘ 
Elimination of the dominant component 
vere . . . . 7 s . 

10. A practical disadvantage of ‘partial transfer’ is that (20), normally, 
entail far more labour than (14). That objection does not apply to a pro- 
cedure which Fox and Green employed in 1941 (2) in combination with 
another device which is explained below. Starting from any guessed mode 

w, (say), they derived from this another mode w, by solving—in conjunc- 
tion with the edge-conditions—the equation 
UW \[I—(Ap) 4 O Jwy, (22) 
in which A as before denotes some convenient multiplier and (Ap) , is the 
al Rayleigh estimate’ deduced from (11) and (12) with w replaced by w,. An 
mee irgument on the lines of section 7 shows that, when 
> [Awe], 
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then, as deduced from (22), 


wp = A > [Ag wet — (Ap) a/Ax}] 


w 


i—Ar)svy When dw, = A#P'w,. (23) 
Thus w, can be regarded as the difference between w, and the result of one 
intensification. (This significance was not remarked in (2).) 

Here we shall term the procedure ‘elimination of the dominant compo- 
nent’. That it achieves this end is apparent from the first expression for 
wy, in (23); for if A,w, (say) is the predominant component of w,, then 


8 


(by “Rayleigh’s principle’) 1—(Ap) ,/A, = 9. 


‘Optimal synthesis’ 

11. In (2) it was used in conjunction with a device which we now 
explain. (It should be emphasized that such combination is not obligatory. 
Either method can be employed without the other.) 

W,, Wy, now denoting any two modes which satisfy the edge-conditions, 
if 


w Wrtyw, (y real) 
the ‘Rayleigh estimate’ derived from w is 


a+2yb+-y?c N 
Ap = - ~ a (say), 
OR d ) A 


where : 
’ 9g 9 “ 
w, dw, dady 


2 


b {| w, Pw, dady {| Uy, Ow, dxdy (24) 


io) 


Uy, Puy, dady 
d {| w, Ow, dady 
f [fee uy, dady I! Uy, Ww, dady 
¢ {| Wy Pw, dxady 
From this expression it follows that Ap is stationary both when y? = « 


and when 


(bd —af) + (cd—ag)y-+ (cf—bg)*)? = 0, (25) 
and that in the latter event it is given by 


Ap = (a+by)/(d+fy) = (b+cy)/(f+ gy). (26) 


When y = 0, 
Ap = a/d = (Ap)4, say, and 40Ap/ey = (Ap) 4{b—f(Ap) ,}/a, 
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so unless b = f(A,), the gradient of A, will not be zero; consequently (if 
(Ap) 4 > ©) some real value of y will give to Ap a positive value nearer than 


/ 


(Ap), to Ay. If (Ap), > 0 it entails no restriction to assume that either 


(Ap) yp > Ar)4 or (Ap)p < A_, (cf. section 6). In either event a closer estimate 


i 


of A, than either (Ap), or (Ag), will be obtained by ‘optimal synthesis’. 


12. Both roots of (25) are real when (as here: cf. section 3) N in (24) is 


necessarily positive ; for on that understanding 


a > 0, c > 0, ac > b’, 
therefore a*c? > ab’c, (i) 
sao d/a—g/c)? d/a+-¢/c)?—4dg/ac 
t\f/b(d/a+¢/c—f/b)—dg/ac} (ii) 
because d/a+¢/c—2f/b)? > 0; 


and on multiplying corresponding sides of (i) and (ii) we have 
(cd—ag)? 1'/bf(cd--ag)—f?ac—b*d¢' 1(bd — af)(cf—bg). 
The modes derived from (25) are conjugate, viz. 


uw. 


w, (say) w Vy Up. Wry (say) w B: 


cvT71 “eB 47 Y2 


for the conditions of such conjugacy are (cf. (9), section 4) 
() || wy, Fwy, dady a b(y, -Yo)+C 
sit according to (24), 


0 || wporwy dxdy d f(y, 


und ieccording to (25) 


(cf—bg)(y,+y2 ag—cd, (cf—bg)y, y. = bd—af. (27) 


Rayleigh’s method in conjunction with intensification 

13. G. Temple (10) and Temple and Bickley (11), also J. J. Koch 
12) in effect, have used intensification in conjunction with Rayleigh’s 
method. By this procedure, when the eigenvalues are all positive, the accuracy 
of his estimate is improved: for as derived from w it is given by (12), and as 
derived from w, it is given by 


Ay )y wl aa say, 
where l | | w, dw, dxdy A? > (AZ/AZ) |. (11)* (28) 
I ( ( wy UW 11 dudy A3 > (A? A?) 


Hence (section 6*) 


(Ar) /Ar | when every A is positive but not (necessarily) otherwise. 
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An alternative to the ‘Rayleigh estimate’ 
14. But all of Rayleigh’s arguments apply to 
pz (say) AST" y (Az) > (A? dz), (13)* (29) 


_ 


whether all of the eigenvalues have one. sign or no, because every 2 > 0. 
That is to say (section 7*) 


pu > the smaller of AZ and A2.,, 
(uz); < we, when (2), is derived from w,. 
Finally, by (12) and (29), 
p2|(Ap)® = 12/11" = [ S (ARAB) P/ (AB). (AD 


= 1 — } [AP ASA, *—AZ1)*]/ > (AZ/AR). > (AR) (30) 


< 1, whatever the signs of the eigenvalues. 


Thus p2 (which of course entails ‘intensification’) has advantages over 
the Rayleigh estimate A;,—more especially in relation to systems like our 
test example. 


15. ‘Optimal synthesis’ (section 11) can be performed on the basis of 


(29). That expression gives for w = w,-+-ywy (y real) 


a+2yb+y"c 
D12)F+°G 


where a, b, c have their expressions (24) and where 


(us/A?),, 


ND = [ [ w', du’, dady \, (24)A 
MEF = | [ w'y dw’, dady 
NG = | [ wy, dw’, dxdy 


w', and w?, denoting functions obtained by intensification of w, and wy on 
the basis of equations like (14). 

Every consequence drawn from (24) has a parallel which can be 
drawn from (24)A, except that (cf. section 10* of (4)) the first but not the 
second of (9) is satisfied by the modes which make ,2 stationary, whereas 


both of (9) are satisfied (cf. section 12) by the modes which make Ap, 
stationary. 
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II. Review or EARLIER WorK ON THE TEST EXAMPLE 


16. In the earlier attack on our test example (2), reliance was placed on a 
use of ‘Rayleigh’s method’ in combination with ‘optimal synthesis’ as 
described in section 11, There, no restriction is imposed on w, or w, except 
that both must satisfy the edge-conditions: in (2), ‘point-relaxation’ was 


employed to relate w, with w, by (22)—no attempt being made to solve 














Fic. 4. 


9° 


(22) exactly, since the aim was merely to obtain a second component 
differing » idely from w,. The objective was strictly practical—to com- 
pute A, in a case where the gravest mode could not be guessed. ‘Intensifica- 
tion’ was not utilized, or ‘partial transfer’; and the alternative (42) to 
Rayleigh’s estimate (Ap) had not at that time been recognized. 

To make the example a drastic test of method, the initial estimate of the 
gravest mode was deliberately made wide of the mark (viz. symmetrical 
both in a2 and y, and of one sign throughout). The first and final nets had 
mesh-side d/4 and d/8 (Fig. 2). On the first, by 12 optimal syntheses the 
mode was greatly altered and the Rayleigh estimate was brought down 
toa nearly stationary value. Fig. 4 (a rough plot of hitherto unpublished 
computations) shows the results of this stage of the work. 

Further syntheses on the finer net brought Ap, down from 7-904 to 5-949, 
a value ‘believed to be correct within about 1 per cent.’ Less exactitude 
was attempted in the determination of the mode (Fig. 12 of (2)—here 
reproduced as Fig. 2); and no proof was adduced to show that this 
approximates to w, (not to some other characteristic mode). 


17. Other theoretical questions were left unanswered: e.g., could the 
accuracy of Fig. 2 have been improved without excessive computational 
labour? is ‘C’, section 1, a feature of this system; and if so, what mode 
is associated with A,? This last question has some practical interest: 
first, since if A, = A, the figure given above (5-949) overestimates 4,; 


secondly, because the mode conditions the best arrangement of ‘stiffeners’. 


5092.24 


Hh 
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These questions prompted an investigation which has resulted both in 
this paper and in (4). At the outset it was evident that any or all of *A’, 
‘B’, and ‘C’, section’ 1, might be features of the test example; that ‘A’ 
would invalidate use of ‘Rayleigh’s method’ except in conjunction with 
some means to an accurate estimation of w,; and that as a consequence 
of ‘B’ either sign might be taken by a ‘Rayleigh estimate’. At first this 
was taken as an indication that Rayleigh’s method must be discarded; but 
subsequently our recognition of an ‘excluded range’ (section 6) showed 
that every positive estimate > A, and so made it likely that the estimate 
5-949 is close, since A, became nearly stationary in the later syntheses 
performed by Fox and Green. 

To this extent their Fig. 4 gained plausibility, despite the italicized 
proviso of section 7. It did not seem possible to prove that w, was a com- 
ponent in their initial guess—and only so, in theory, would it contribute 
to the modes which they derived from it by a method shown (section 10) 
in effect to be ‘intensification’; but the likelihood of such failure is reduced 
by the fact that (22) was not solved exactly. None the less it was felt that 
other methods should be tried, to attain their results, if possible, by a 
different route. 


[1]. New Meruops APPLIED TO THE TEST EXAMPLE 
(1) ‘Exploration’ with an applied point-loading 
18. This method rests on the theoretical expression for the response (w) 
to a concentrated force applied at some chosen point P when the intensity 
of the destabilizing forces (edge-thrusts) is specified by A. In (16), section 8, 
let Z, be zero except within an elemental area at P, and let 


| | Z, dady for the elemental area = 1. 
Then, on multiplying (16) by w,, and integrating over the whole area of the 
plate, we obtain 


n 


(w,)p | | w,, Zo dady | | w, Pw, dxdy = A,, simply, 


when its series > [A,w;]| is written for wy. Accorditigly in this instance 
— % . « 


(17) become 


Wo z [(w,.)p We], ) (31) 
w > | (wy) p w,/(1—A/A;,) | | 
and show that the response at the ‘load-point’ P is 
(w)p = > [(w,)p/(1—A/A,) 
I > [evn r)| ) (32) 


> [(w;,)2], when A = of 
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19. It follows from (32) that when A = A, (that is, when one of the eigen- 
values is attained) a zero force at P will make (w), unity. The second of 
(32) shows this force to be positive and finite when A = 0; so, on the side 
\ > 0 it first becomes zero when A = A,, on the side A < 0 it first becomes 


| zero when A = A_,. If, then, w can be determined for chosen values of A, 
i 





a ae 


~ TFOR CE at P | 
a |(for 2p 1,000) 
(a)| 0 14,639 


+— + 


b)| 1:6 10,047 











(c)| 3-2 4.783 
(d)| 4-24 278 










\(e)| 4-2 50 | 
~s |} —_+——__+—____, 
| 5 (f) | 4-206 13 
()| 4-208 | 1 | 
| oe aes 4 - —— | 
(e) 
= woh >A 
| 0 1 2 3 (d) 5 
Fic. 5 
by plotting the force required when (w), = 1 (a quantity which we here 
denote by pp) we can derive A, as the smallest positive A for which pp = 0. 


Like ‘partial transfer’ the method effects a bias towards modes and eigenvalues 
of the wanted sign. 

Here too, however, theory imposes a proviso: in order that the method may 
jield wy, (w,)p must not be zero. Thus here too, in theory, ‘A’ may be a 


uuse of failure; but the chance of this is small—and will be made still 





smaller if ‘exploration’ is repeated with the unit force applied at a different 
point. Whether the method is feasible can of course be settled only by trial. 
It entails use of a different ‘pattern’ at every nodal point; but precision is 
needed only in the final stages, when A has a value close to A, and pp = %. 
Application to the test example (Figure 1) 

20. Fig. 5 illustrates the procedure as applied to the coarse-net system. 
Each point represents a relaxational solution in which, with the same 
displacement imposed at P, residuals were brought sensibly to zero 
except at the boundary and at the ‘load-point’ P. On that understanding 
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(cf. section 19) the residual at P is a measure of —pp. Only relative values 
are required to determine A. 
(i) A curve drawn through (a), (6), and (c) suggested a zero force for 
A = 4-24; 
(ii) trial of this value led to solution (d); 
(iii) a new curve drawn through points (a), (b), (c), and (d) suggested a 
zero for A = 4:20; 
(iv) trial of this value led to solution (e). 
Thereafter, linear interpolation could be used: between (d) and (e) to 
suggest a value (4-206) which led to (f), and between (e) and (/) to suggest 
a value (4-208) which led to (g). This last solution entailed a force at P 
which was comparable with the other neglected residuals: consequently 
the mode obtained in (g) could be accepted, and this gave 4-2073 as the 
estimate of A,. 
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21. Fig. 6, which shows the mode as thus deduced, may be compared 
with Fig. 4, which gives the mode obtained by Fox and Green. The general 
similarity of the diagrams, and the agreement of the corresponding Ray- 
leigh estimates (4-288 and 4-2073), suggest that both relate to the gravest 
mode; but the net is too coarse to define this with precision, or to justify 
replacement of the integrals by summations. (Fox, using standard formulae 
of integration, found on halving the mesh-side that A rose at once from 
4-288 to 7-904.) 


Extension of the method 

22. Inthe unlikely event (section 19) of P’s being nodal for w,, A, will not 
be altered by the imposition of a constraint at P, so an exploration made 
with P restrained and another point Q loaded will yield the same eigenvalue 
as before; but if (w,)p 4 0 then, by a theorem due to Rayleigh ((5), section 


88), a constraint at P will cause A, to rise. We may exclude as almost 
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3 unthinkable the possibility that both P and Q are nodal for w,; and on 

that understanding two explorations (one as above; the other with P loaded 

r and Q restrained) will indicate whether ‘C’, section 1, is a feature of the 
system. 

For suppose that A, = A, exactly. Then two distinct modes are associated 

with A,, and in some combination of them Q is nodal; therefore if (as in 


section 19) pp denotes the force required to make (w)p unity, and if pq 
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stands for the corresponding force of restraint at Q, gp as well as pp will 
be zero when A = A, as obtained before. If A, > A, but only by a small 
quantity, the constraint at Q (which we have assumed not to be nodal for 
d w,) will cause pp to vanish for a slightly higher A, and gp will not vanish 
u simultaneously. 
i 
st Application to the test example 
'V 23. Figs. 7 present the results of a double exploration on the same net 
Le (a = })as before. In (a), with Q(Fig. 5) restrained a fixed displacement was 
m imposed at P, and pp = © when A = 4-27. Comparing this value with 
| 4-208 as deduced from Fig. 5 (ef. section 20), we conclude that (for this net) 
(i) A, 4 A,—otherwise the estimate of A, would not have been raised 
(cf. section 22); 
‘ (ii) P is not nodal in w,—otherwise (unless Q too is nodal in that mode) 
le Figs. 5 and 7(a) would have yielded widely different estimates of A,; 
- (iii) Q is not nodal in w, (since the estimate of A, is raised as a consequence 
yn 


of the restraint at Q); but lies close to a nodal line (since the rise of 


A, is less than 2 per cent.). 
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Both (ii) and (iii) are confirmed by Fig. 6 if (as in section 21) that diagram 

is assumed to approximate to w. 
In (6), with point P restrained a fixed displacement is imposed at Q, and 

qq == 0 when A = 4-44. This increased rise gives further support to (iii). 
24. In theory a ‘Maxwell relation’ should be satisfied—namely 


pq Gp (33) 
and Fig. 5 should be deducible from Figs. 7 by the formula 
pp in Fig. 5 = (pp—pq?/qq) as computed from Figs. 7. (34) 
(33) is confirmed less exactly than was expected; but (34) is very closely 
realized when the product q.qp is taken in place of pq?. Both pq and gp 
come to zero before either pp or 7. 
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That the curves of pp and qq in Figs. 5 and 7 are nearly linear is evidence 


that w, predominates in (32). It does so because P and Q lie (Fig. 5) in 


the region where the destabilizing forces are large: if either or both had 


‘ 
<4 


lain near the top right-hand corner in Fig. 1, then initially ‘negative’ modes 
would have predominated and pp and gq would have increased with Ad, 
falling only (and then rapidly) to zero as A > A,- Results of this kind were 
obtained in an earlier exploration, when P had the same position as 
before but Q was on the centre-line: then (cf. Figs. 8), though jp still fell 
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m fairly steadily, the variation of gq was quite different. (The “Maxwell 


relation’ (33) is closely satisfied.) 











( i 
“ Practical aspects 

25. So far as it can be tested by one example, ‘exploration’ seems to meet 

all practical requirements when (as on the coarser net) use of a different 
>) relaxation pattern’ at every point is feasible. Each solution yields a 
good starting assumption for the next, and (if the point of loading is chosen 
1) suitably) no trouble comes from ‘negative’ components. But (section 21) 
a a finer net is needed to define the mode and thereby to evaluate A, with 
3 precision; and no method can be recommended which entails a different 
p pattern at each of the 273 ‘relaxable’ nodes in Fig. 2. Thus ‘exploration’ 
has value only in the initial computations on the coarser net: there, it 
provides an estimate of the gravest mode which has no appreciable 
negative’ component and accordingly, if made the starting assumption 
| on the fine net, will almost certainly obviate trouble due to ‘A’ and/or 

| ‘B’, section 1. 

But nothing found on the coarse net can establish that two or more of 
the gravest eigenvalues will not be nearly or exactly equal in the fine-net 
system. If A, = A, exactly, any combination of w, and w, may be taken 
as the gravest mode and will not be altered by intensification; but if A, is 
slightly greater than w,, then in theory intensification will reduce the 
component w,, but very slowly, so in practice it may fail to yield w, on 
account of inescapable computational errors. 

26. On the fine net Fox and Green used optimal synthesis in conjunction 
with the method of section 12 (there shown to be equivalent to ‘intensifica- 
tion’) and ‘after point liquidation had been carried as far as was deemed 
practicable’ ((2), section 32). Probably they would have found these 
methods insufficient if they had sought to determine the mode as accurately 
as the eigenvalue; for while one intensification of their accepted mode 
(Fig. 2) reduced A,;, by only 0-1 per cent.,+ both it and subsequent intensifica- 
tions entailed appreciable changes in the mode. Other attempts based on 
various starting assumptions agreed in suggesting that A, lies close to 5-94; 

' but in no instance was a quite definite mode obtained, and as a rule each 

: intensification (in the latter stages) resulted in a small but appreciable and 

| nearly invariant alteration (in the notation of section 7, A, w,—w and 

: Aly = Wy; — wv Were small and nearly identical). Moreover, occasionally A, 

4 appeared to rise 

" 

: Name from 5:9464 to 5-9413. The value (5-949) which Fox and Green recorded is 
believed to relate to a mode slightly different from Fig. 2; but their conclusion that it is 


| orrect within about | per cent Ss ipported by these figures. 
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It thus appeared that ‘intensification’—the only method yet available 
on the finer net 





might fail as applied to systems which have the features 
‘A’ and *B’ and may have the feature ‘C’, section 1. (The fact—shown by 
Fig. 7—that A, and A, are not equal for the coarse-net system cannot be 
adduced as an argument against their being equal when the number of the 
eigenvalues is much larger.) 


(2) Introduction of an ‘ancillary restraint’ 

27. The possibility that A, = A, exactly may be tested by introducing an 
‘ancillary restraint’. Rayleigh ((5), section 92a) has shown that application 
of a constraint (unless at a point which is nodal: there it has no effect) 

(i) reduces by one the number of the eigenvalues ; 

(ii) raises each eigenvalue to a new value intermediate between its old 

value and that of the next higher eigenvalue. 
If originally the two gravest eigenvalues were exactly equal, this means 
that one is conserved, but only one: the second is made greater than before. 
If more than two are equal, by employing more than one constraint we can 
again make the gravest eigenvalue unique; and for the modified system it 
appears (from one trial) that intensification is effective. 

Assuming an ‘ancillary restraint’ to prevent displacement of a point P 
in the region where the destabilizing forces are large, by intensification we 
can find the gravest mode and eigenvalue (w, and A,, say) of the system 
as thus modified, also the force exerted by the restraint. If this is zero, then 
the restraint is inoperative, and A,, w, may be identified with the wanted 
A,, w,: if it is large, then the restraint entails a large increase of stability, 
therefore (A~—A,) and a fortiori (A,—A,) must be considerable (for we have 
seen that A, lies between A, and A,). The latter result is unlikely, since 
unless A, = A, intensification should not have failed to yield w,. More 
probably, what will result are a mode (w,) entailing a small force at P and 
an eigenvalue (A,) slightly greater than A,. 


Theoretical aspects 

28. G. V., applying this procedure to the test example, took as her 
starting assumption a combination, making w, zero, of two modes found 
earlier and believed to be fair approximations to w,. From this she was 
able by intensification to derive an acceptable solution wo and a value 
6-0032 for Ap. 

According to (31) of section 18, 


Woe = > [ (wz) p w,/(1—Ag/A,)], Wy = > [(wx)p wx] | 
so Wo (say) = We—Wy = Ao DY [ (Wy) p Wel (A—Ac)] 


, (35) 
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Each of we and w, has a singularity due to unit forces at P, so W, has no 
singularity; and the first terms dominate its series (35) by reason of the 
| smallness of (A,—Ag) and (A,—Ac) compared with other denominators 
always provided that (w,)p, (w,)p are not zero. This is made reasonably 
certain by P’s location in a region where the destabilizing forces are large. 
In G. V.’s test, P had the position shown in Fig. 5.) 
Two intensifications of W,, in accordance with 
OW, \ o'W,, oW, = Ad'W, 
yi ld modes 


WW 
VW 


( \ ~ > | (W;.) p Wi, Aj(Aj, A-)] | 


: 7 (36) 
o = Bg ¥ [(we)p wel ARA—Ac)] J 


which also have no singularity at P and in which w,, w, are still more pre- 





dominant. The Rayleigh estimate for W, < Ag, the Rayleigh estimate for 
| w,: but will not be much smaller if the force of restraint is small. 

G. V.’s computations gave a value slightly less than 6-002. This being 
larger than some earlier estimates, it was to be expected that W., W. would 
differ sensibly—and this proved to be the fact. They had the characteristic 
noticed in section 26: namely, that A,, A,,; as there defined were small and 
nearly identical. The suspected reason was near-equality of A, and A,; and 
the final problem, in the test example, was seen to be closer estimation 


ol the oravest mode. 


The final problem: estimation of the gravest mode 

29. At this point our collaboration was interrupted and it was left to 
the senior author to make a final assessment of results. That task led him 
to some theorems which were presented in (4) and, in respect of the 


test example (sections 17—19*), to the conclusion that e7ther A, is restricted 


t 


V the doubk inequality 


5-9648, <A 5-9967, (38*) (37)A 


1 


(and in that event the earlier estimate, 5-949, of Fox and Green (section 15) 
isin error by not more than 0-8 per cent.) or A,, which is very nearly equal 
to A, lies within the range 
596480. <A < 6-03943, (39)*, (37) B 
and then A, may be less than the lower of those limits. In this second event 
the estimate of Fox and Green is likely to be still closer. 
The inequalities (37), A and B, are based on values which G. V., identi- 
fying w, w,, w;; with Wo, Wo, Wo, and giving to A the value (6-002) which 
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she had found for A, (section 28), computed for the integrals 
I = || wow dedy = ¥ [Az] 
I’ [ ( w dw, dady = A > | AZ/A,] 
1" = [| wow dedy = A* Y [ABA] = [[ wy Pw, dady + —(15)* (38) 
I” = || wy Swy dedy = A ¥ [ARR] 
[= ff wy, Sw, dedy = AS [ARIAL] 


(Her identification of A with Ap entailed exact equality of J and I’ (cf. (12), 
section 6), and her other estimates were 
I" /I 1-0012,, iad 1-0019., T'v/T = 1-0029, (32*) (39) 
a dropped figure such as 3 denoting some value between 2 and 4 but 
possibly nearer to 2 or 4 than to 3.¢) In addition to these restrictions on d,, 
(4) drew some conclusions in regard to modes: e.g. (section 18*) that as 
computed from (12) for 


ww, wtyw, (= WetyWe) (23)* (40) 

A,/A has a maximum value — 2-464 15039 when y —0°99914 and | 
(41) 

a minimum value 0-999 140032, when y = 2-47 


and that as computed from (30) for w, 
ue/A* has a maximum value (2-081 307 18)? when y 0-99907 and | 
a minimum value (0-999 11809)? when y = 3-40. 


(33)* (42) 


30. Only the second of (42) was utilized, (4) being only concerned with 
eigenvalues; but all of (41) and (42) serve to elucidate the final problem 
(section 29) when taken in conjunction with (37) A. ThusApand pz as com- 
puted for w, take maximum values when y - | with an error of the 
order of 0-1 per cent.: i.e. when to this order of approximation the mode is 

wW—Wy A,, as defined in section 26, 
> [A;,w,(1—A/A,)] when w= } (A;,w;,). (48) 
Now in (4) (37) was obtained as 
0-99380, < A,/A < 0-999118 (37)* (44) 


+ Strictly, the computations in (4) should have been extended to allow for this possibi ity. 
But thereby the paper would have been greatly lengthened, and it is unlikely that the 


limits would have been much altered. 
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with A = 6-002 (section 29). This shows that w, makes a negligible con- 
tribution to A,, therefore A? has a negligible influence on the values which 
A, and 2 take for A;; and the arguments of section 7 and section 14 then 
| serve to deduce from (41) and (42) that 





A_,|/A 2-46415, A,/A < 2-08131. (45) 
31. Wecan go further if in each of the series (43) and (38) three terms, at 
most, are taken as significant, so that we may write 
w= A,w,+A,u,+A,u,, I = Aj+A?-+ A...., ete. (46) 
For then elimination of A? between successive pairs of (38) yields four 
relations between J, A?, and A?>—namely, 
' (A, xB, 2?C, x#°D)/ (1,7, r?, r3)(1—r)A2-+ (1, 8, s?, s)(1—s)A? 
where ; 
A A, \ / A, A, s A, A, (47) 
A | e (Say), B e—ra 
LS 
C e+a—ab, D «e+b—zc 
) ind when | 
lta I" /I l--b ‘al @ lc iv/T }, (48) 
so that 104 x (a, b, c) 12°,, 19°,, 29°53 
according to (39); and all of (47) are closely satisfied when 
} A, A s 0-4a | 
r)A®/T l—a(1-+a/1-4) = 1—z/0-9992 }. (49) 
12/] a/(1-4)? 6-2 x 10-4 
since 7 l (by hy pothesis ind x ( A, \) is restricted by (44), (49) give 
n ; — ‘ t a y - 
real values both to A, and to A,. A? is very small in relation to J (and to A?), 
n , : ; , 
| but A? can havea value which is comparable, since (1 —r) is small like (1—2). 
1e nn - ° ° 
Thus the values (39) afford indications that 
IS 
(a2) another characteristic mode (w,) is associated with an eigenvalue (A,) 
which is very little greater than A,. 
° b) the mode (W,,) from which (39) were derived contains an appreciable 
‘negative’ component, associated with an eigenvalue (A_,) nearly 
{ equal to A, (0-4) 2-5A,. The amplitude, A,, of this negative 
component of W,, is small (of the order 10~-*) in relation to A, and A,; 
. but in A, its amplitude is sensible because the other components are 
reduced almost to zero. 
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Undue weight must not be given to the indications, since they are based 
on a merely plausible hypothesis and on computed values which have 
uncertain accuracy. But the indications are consistent with (45) and with 
the observations recorded in section 28, also with the second of the possibili- 
ties that were contemplated in section 29; and they give helpful guidance 
towards a solution of the final problem (section 29). (For example, given 
a close estimate of A, (= A,/s), the negative component of W;, can be elimi- 
nated by the method of section 10.) 
32. The indications, moreover, are unambiguous—to the extent that 
the values (39) are accurate. The first two of (47) require that 
(r—s)(1—r)A? = (xB —s8A)) 
ds (50) 

(s—r)(1—s)A?2 (xB—rA)I} 

and the last two, when these expressions are substituted, give 

(r+s)¢aB—rsA = 2?C 
(r+s)?2B—rs{~aB+(r+s)A} = 23D, 

whence (B?—AC)(r+-s) = 2(BC—AD) 

(B?—AC)rs = 2?(C?— BD) 


(51) 


Hence r and s can be computed. 
The inequality (44) implies that 
8 < «x 104 < 62, 
so all of «, a, b, c are of order 10-3 or less, and terms of the third and higher 
order in those quantities are negligible in comparison with a?, ab.,..., ete. 
On that understanding, when A, B, C, D are given their expressions (47), 
B?— AC = a?—«(3a—b) | 
BC— AD = a(b—a)+«(c—2b) 

2—BD = b*—a(b-+c—a)+ «(c+ 3a—3b) | 
with sufficient approximation. These quantities vanish for values of € x 10* 
which differ, but which all lie between 8-5 and 9-5. In that narrow range 
as computed from (51) and (52), (r+-s) and rs vary rapidly—and take values 
excluded by the requirements |r| < 1, |s| < 1: but all of (B?—AQC)...., etc., 
vanish for the same value of e+ when changes are made in a, b, and c which 


+ The necessary and sufficient condition is 
a(c— 2a) = (b—a)?, 
and is closely satisfied when, for example, (48) are replaced by 
10* x (a, b, c) = 12-2, 19-5, 28-8. 
Then 10! x € 8-8, approximately, according to all of (52). 
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hardly go outside their margins of uncertainty (cf. section 29); and (52) then 

, indicate that 

(7--8)/x 0-6, rs/x* = _— 

therefore rit, 1% = i, 0-4, 


i.e. A= A,/z, A, = —25A,/x 


s 


(53) 


These are the values recorded in (49), and when A? and A? have the values 
given there the first two of (47) are satisfied exactly, the last two very 
nearly if « lies within the limits imposed by (44). The last two of (53) are 








compatible with (45); and the first of them indicates very near equality 
of A, and A, (thus explaining the near-identity of A; and A,,t), the second 
explains the negative maximum of Ap which is given in the first of (41). 
Conclusion 
33. Finally—and consistently with the indications of section 31—the 
| second of (42) shows that a better approximation to w, than W, is the mode 
w = (W.+3-4W,,)/4-4 (54) 
of which contours are given in Fig. 9. This differs sensibly from the mode 
obtained by Fox and Green (Fig. 2), nodal lines of which are reproduced 
in Fig. 9: yet both modes are associated with eigenvalues that differ from 
5-97 by less than | per cent. 
Thus although the ‘final problem’ has not been solved completely, what 
has been done is of nearly equal value from a practical standpoint. For 
1 Figs. 2 and 9 are different combinations of w, and w, (each containing, 
probably, a small admixture of w_, and of other unwanted modes); and all 
combinations are almost equally likely to occur (since A, and A, have been 
shown to be nearly equal), so the conclusion in regard to stiffener (section 17) 
is evident: the region of maximum waving has been delimited, but the 
) phase’ of the waving is nearly indeterminate. 
34. All three of the difficulties ‘A’, “B’, and ‘C’, section 1, have been found 
yA 


to characterize the test example; and all for practical purposes have been 
surmounted, by methods which should have very wide application. We 
have shown, in sections 18—24 and 27, how to ensure the presence of w, 
ina guess’ which is to be improved by intensification: in (4), how, at the 
cost of two intensifications, to impose both upper and lower limits on d,: 
in sections 31—32, how from two intensifications to draw helpful indications 
regarding the three modes w,, w., w 


Ll? yy 


The cofactor of A, w, is —7 in Ay, —n(1—7) in Ay, where 7 = 1—A/A,; and these 


quantities are nearly equal when A, A, /2 A. 
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\ word should be added concerning our reference to methods by descrip- 

tive titles instead of names. We are in no doubt about the attribution of 

| ‘Rayleigh’s principle’, but except in this instance we have hesitated to 
pronounce on questions of priority. It would seem that some of the 
methods we have described were devised by two or more authors indepen- 
dently—as a rule, for the solution of particular problems; and that their 
general applicability was recognized only later, by expositors. We have 
found attribution specially difficult in respect of the many methods—e.g., 
Holzer’s’—which in effect satisfy the equation for normal modes at all 
but one point of the system considered. In our view the purpose of this 
treatment is material. It was proposed, for example, by Rayleigh in the 
second edition (1926) of his Theory of Sound (section 88) as a means of 
approximating to eigenvalues, and we think it unlikely (since he gave no 


reference) that he derived his notion from Holzer s book (13: 1921): 








| certainly this did not influence Bairstow and Stedman (14), who proposed 


the method in 1914 (for struts 


This paper in effect employs the treatment, but with emphasis on the 
position of the loaded point, since its concern is rather with the mode than 
with the eigenvalue. We hesitate, in this matter as elsewhere, to make 


claim to any priority: all we would say is that sections 8-12, 14, 18-24, 


>.) 


and 32, also (4), contain matter which to the best of our belief is new. 


35. Grateful acknowledgement is due to the Clothworkers’ Company, 


for a grant to Imperial College which made Miss Vaisey’s collaboration 


possible: also to Miss M. Evenett, of that College, for assistance with the 
diagrams and in computations, requiring a desk machine, which have been 
needed since that collaboration ended (section 29). (Neither she nor 
Miss Vaisey, in the nature of the case, has responsibility for the conclusions 


ot sec tions 29-33 
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